Image for :meth`.Surface.evaluate_cartesian_multi` docstring.
508508
def surface_evaluate_cartesian_multi(surface, points):
    """"""
    if NO_IMAGES:
        return

    ax = surface.plot(256)
    ax.plot(
        points[0, :], points[1, :], color="black", linestyle="None", marker="o"
    )
    delta = 1.0 / 32.0
    font_size = 18
    ax.text(
        points[0, 0],
        points[1, 0],
        r"$w_0$",
        fontsize=font_size,
        verticalalignment="top",
        horizontalalignment="right",
    )
    ax.text(
        points[0, 1] + 2 * delta,
        points[1, 1],
        r"$w_1$",
        fontsize=font_size,
        verticalalignment="center",
        horizontalalignment="left",
    )
    ax.text(
        points[0, 2],
        points[1, 2] + delta,
        r"$w_2$",
        fontsize=font_size,
        verticalalignment="bottom",
        horizontalalignment="left",
    )
    ax.axis("scaled")
    ax.set_xlim(-3.125, 2.375)
    ax.set_ylim(-0.25, 2.125)
    save_image(ax.figure, "surface_evaluate_cartesian_multi.png")
881448
def curve_evaluate(curve):
    """"""
    if NO_IMAGES:
        return

    ax = curve.plot(256)
    points = curve.evaluate_multi(np.asfortranarray([0.75]))
    ax.plot(
        points[0, :], points[1, :], color="black", linestyle="None", marker="o"
    )
    ax.axis("scaled")
    ax.set_xlim(-0.125, 1.125)
    ax.set_ylim(-0.0625, 0.5625)
    save_image(ax.figure, "curve_evaluate.png")
368357
def mosaic_inline(self, imagelist, bg_ref=None, trim_px=None,
                      merge=False, allow_expand=True, expand_pad_deg=0.01,
                      max_expand_pct=None,
                      update_minmax=True, suppress_callback=False):
        
        # Get our own (mosaic) rotation and scale
        header = self.get_header()
        ((xrot_ref, yrot_ref),
         (cdelt1_ref, cdelt2_ref)) = wcs.get_xy_rotation_and_scale(header)

        scale_x, scale_y = math.fabs(cdelt1_ref), math.fabs(cdelt2_ref)

        # drop each image in the right place in the new data array
        mydata = self._get_data()

        count = 1
        res = []
        for image in imagelist:
            name = image.get('name', 'image%d' % (count))
            count += 1

            data_np = image._get_data()
            if 0 in data_np.shape:
                self.logger.info("Skipping image with zero length axis")
                continue

            # Calculate sky position at the center of the piece
            ctr_x, ctr_y = trcalc.get_center(data_np)
            ra, dec = image.pixtoradec(ctr_x, ctr_y)

            # User specified a trim?  If so, trim edge pixels from each
            # side of the array
            ht, wd = data_np.shape[:2]
            if trim_px:
                xlo, xhi = trim_px, wd - trim_px
                ylo, yhi = trim_px, ht - trim_px
                data_np = data_np[ylo:yhi, xlo:xhi, ...]
                ht, wd = data_np.shape[:2]

            # If caller asked us to match background of pieces then
            # get the median of this piece
            if bg_ref is not None:
                bg = iqcalc.get_median(data_np)
                bg_inc = bg_ref - bg
                data_np = data_np + bg_inc

            # Determine max/min to update our values
            if update_minmax:
                maxval = np.nanmax(data_np)
                minval = np.nanmin(data_np)
                self.maxval = max(self.maxval, maxval)
                self.minval = min(self.minval, minval)

            # Get rotation and scale of piece
            header = image.get_header()
            ((xrot, yrot),
             (cdelt1, cdelt2)) = wcs.get_xy_rotation_and_scale(header)
            self.logger.debug("image(%s) xrot=%f yrot=%f cdelt1=%f "
                              "cdelt2=%f" % (name, xrot, yrot, cdelt1, cdelt2))

            # scale if necessary
            # TODO: combine with rotation?
            if (not np.isclose(math.fabs(cdelt1), scale_x) or
                not np.isclose(math.fabs(cdelt2), scale_y)):
                nscale_x = math.fabs(cdelt1) / scale_x
                nscale_y = math.fabs(cdelt2) / scale_y
                self.logger.debug("scaling piece by x(%f), y(%f)" % (
                    nscale_x, nscale_y))
                data_np, (ascale_x, ascale_y) = trcalc.get_scaled_cutout_basic(
                    data_np, 0, 0, wd - 1, ht - 1, nscale_x, nscale_y,
                    logger=self.logger)

            # Rotate piece into our orientation, according to wcs
            rot_dx, rot_dy = xrot - xrot_ref, yrot - yrot_ref

            flip_x = False
            flip_y = False

            # Optomization for 180 rotations
            if (np.isclose(math.fabs(rot_dx), 180.0) or
                np.isclose(math.fabs(rot_dy), 180.0)):
                rotdata = trcalc.transform(data_np,
                                           flip_x=True, flip_y=True)
                rot_dx = 0.0
                rot_dy = 0.0
            else:
                rotdata = data_np

            # Finish with any necessary rotation of piece
            if not np.isclose(rot_dy, 0.0):
                rot_deg = rot_dy
                self.logger.debug("rotating %s by %f deg" % (name, rot_deg))
                rotdata = trcalc.rotate(rotdata, rot_deg,
                                        #rotctr_x=ctr_x, rotctr_y=ctr_y
                                        logger=self.logger)

            # Flip X due to negative CDELT1
            if np.sign(cdelt1) != np.sign(cdelt1_ref):
                flip_x = True

            # Flip Y due to negative CDELT2
            if np.sign(cdelt2) != np.sign(cdelt2_ref):
                flip_y = True

            if flip_x or flip_y:
                rotdata = trcalc.transform(rotdata,
                                           flip_x=flip_x, flip_y=flip_y)

            # Get size and data of new image
            ht, wd = rotdata.shape[:2]
            ctr_x, ctr_y = trcalc.get_center(rotdata)

            # Find location of image piece (center) in our array
            x0, y0 = self.radectopix(ra, dec)

            # Merge piece as closely as possible into our array
            # Unfortunately we lose a little precision rounding to the
            # nearest pixel--can't be helped with this approach
            x0, y0 = int(np.round(x0)), int(np.round(y0))
            self.logger.debug("Fitting image '%s' into mosaic at %d,%d" % (
                name, x0, y0))

            # This is for useful debugging info only
            my_ctr_x, my_ctr_y = trcalc.get_center(mydata)
            off_x, off_y = x0 - my_ctr_x, y0 - my_ctr_y
            self.logger.debug("centering offsets: %d,%d" % (off_x, off_y))

            # Sanity check piece placement
            xlo, xhi = x0 - ctr_x, x0 + wd - ctr_x
            ylo, yhi = y0 - ctr_y, y0 + ht - ctr_y
            assert (xhi - xlo == wd), \
                Exception("Width differential %d != %d" % (xhi - xlo, wd))
            assert (yhi - ylo == ht), \
                Exception("Height differential %d != %d" % (yhi - ylo, ht))

            mywd, myht = self.get_size()
            if xlo < 0 or xhi > mywd or ylo < 0 or yhi > myht:
                if not allow_expand:
                    raise Exception("New piece doesn't fit on image and "
                                    "allow_expand=False")

                # <-- Resize our data array to allow the new image

                # determine amount to pad expansion by
                expand_x = max(int(expand_pad_deg / scale_x), 0)
                expand_y = max(int(expand_pad_deg / scale_y), 0)

                nx1_off, nx2_off = 0, 0
                if xlo < 0:
                    nx1_off = abs(xlo) + expand_x
                if xhi > mywd:
                    nx2_off = (xhi - mywd) + expand_x
                xlo, xhi = xlo + nx1_off, xhi + nx1_off

                ny1_off, ny2_off = 0, 0
                if ylo < 0:
                    ny1_off = abs(ylo) + expand_y
                if yhi > myht:
                    ny2_off = (yhi - myht) + expand_y
                ylo, yhi = ylo + ny1_off, yhi + ny1_off

                new_wd = mywd + nx1_off + nx2_off
                new_ht = myht + ny1_off + ny2_off

                # sanity check on new mosaic size
                old_area = mywd * myht
                new_area = new_wd * new_ht
                expand_pct = new_area / old_area
                if ((max_expand_pct is not None) and
                        (expand_pct > max_expand_pct)):
                    raise Exception("New area exceeds current one by %.2f %%;"
                                    "increase max_expand_pct (%.2f) to allow" %
                                    (expand_pct * 100, max_expand_pct))

                # go for it!
                new_data = np.zeros((new_ht, new_wd))
                # place current data into new data
                new_data[ny1_off:ny1_off + myht, nx1_off:nx1_off + mywd] = \
                    mydata
                self._data = new_data
                mydata = new_data

                if (nx1_off > 0) or (ny1_off > 0):
                    # Adjust our WCS for relocation of the reference pixel
                    crpix1, crpix2 = self.get_keywords_list('CRPIX1', 'CRPIX2')
                    kwds = dict(CRPIX1=crpix1 + nx1_off,
                                CRPIX2=crpix2 + ny1_off)
                    self.update_keywords(kwds)

            # fit image piece into our array
            try:
                if merge:
                    mydata[ylo:yhi, xlo:xhi, ...] += rotdata[0:ht, 0:wd, ...]
                else:
                    idx = (mydata[ylo:yhi, xlo:xhi, ...] == 0.0)
                    mydata[ylo:yhi, xlo:xhi, ...][idx] = \
                        rotdata[0:ht, 0:wd, ...][idx]

            except Exception as e:
                self.logger.error("Error fitting tile: %s" % (str(e)))
                raise

            res.append((xlo, ylo, xhi, yhi))

        # TODO: recalculate min and max values
        # Can't use usual techniques because it adds too much time to the
        # mosacing
        #self._set_minmax()

        # Notify watchers that our data has changed
        if not suppress_callback:
            self.make_callback('modified')

        return res
566426
def image(self):
        """"""
        events = self.event_table
        skycoord = SkyCoord(events['GLON'], events['GLAT'], unit='deg', frame='galactic')
        pixcoord = PixCoord.from_sky(skycoord=skycoord, wcs=self.wcs)

        shape = self.config['shape']
        bins = [np.arange(shape[0] + 1), np.arange(shape[1] + 1)]
        sample = np.vstack((pixcoord.y, pixcoord.x)).T
        data, _ = np.histogramdd(sample=sample, bins=bins)
        data = data.astype('float32')

        header = self.wcs.to_header()
        return fits.ImageHDU(data=data, header=header, name='image')
179152
def _spherical_to_cartesian(cls, coord, center):
        
        r, theta, phi, r_dot, theta_dot, phi_dot = coord
        x = r * cos(phi) * cos(theta)
        y = r * cos(phi) * sin(theta)
        z = r * sin(phi)

        vx = r_dot * x / r - y * theta_dot - z * phi_dot * cos(theta)
        vy = r_dot * y / r + x * theta_dot - z * phi_dot * sin(theta)
        vz = r_dot * z / r + r * phi_dot * cos(phi)

        return np.array([x, y, z, vx, vy, vz], dtype=float)
1312973
def bbox(img):
    
    rows = np.any(img, axis=1)
    cols = np.any(img, axis=0)
    rmin, rmax = np.where(rows)[0][[0, -1]]
    cmin, cmax = np.where(cols)[0][[0, -1]]

    return rmin, rmax, cmin, cmax
362556
def correct(self, img):
        '''
         
        '''
        print("CORRECT PERSPECTIVE ...")
        self.img = imread(img)

        if not self._homography_is_fixed:
            self._homography = None
        h = self.homography

        if self.opts['do_correctIntensity']:
            tf = self.tiltFactor()
            self.img = np.asfarray(self.img)
            if self.img.ndim == 3:
                for col in range(self.img.shape[2]):
                    self.img[..., col] /= tf
            else:
                self.img = self.img / tf
        warped = cv2.warpPerspective(self.img,
                                     h,
                                     self._newBorders[::-1],
                                     flags=cv2.INTER_LANCZOS4,
                                     **self.opts['cv2_opts'])
        return warped
608292
def _get_cartesian_edge_set(self):
        
        # Get projection space for cartesian projection
        edge_sets = numpy.vstack(self.edge_set)
        west, east, north, south = utils.get_spherical_bounding_box(
            edge_sets[:, 0],
            edge_sets[:, 1])
        self.proj = utils.OrthographicProjection(west, east, north, south)

        for edges in self.edge_set:
            # Project edges into cartesian space
            px, py = self.proj(edges[:, 0], edges[:, 1])
            # Store the two end-points of the trace
            self.cartesian_endpoints.append(
                numpy.array([[px[0], py[0], edges[0, 2]],
                             [px[-1], py[-1], edges[-1, 2]]]))
            self.cartesian_edges.append(numpy.column_stack([px, py,
                                                            edges[:, 2]]))
            # Get surface length vector for the trace - easier in cartesian
            lengths = numpy.sqrt((px[:-1] - px[1:]) ** 2. +
                                 (py[:-1] - py[1:]) ** 2.)
            self.length_set.append(lengths)
            # Get cumulative surface length vector
            self.cum_length_set.append(
                numpy.hstack([0., numpy.cumsum(lengths)]))
        return edge_sets
1121693
def img(self):
        ''''''
        SlipThumbnail.img(self)

        if self.rotation:
            # rotate the image
            mat = cv.CreateMat(2, 3, cv.CV_32FC1)
            cv.GetRotationMatrix2D((self.width/2,self.height/2),
                                   -self.rotation, 1.0, mat)
            self._rotated = cv.CloneImage(self._img)
            cv.WarpAffine(self._img, self._rotated, mat)
        else:
            self._rotated = self._img
        return self._rotated
533006
def map(self, coords):
        
        m = np.empty(coords.shape)
        m[:, :3] = (coords[:, :3] * self.scale[np.newaxis, :3] +
                    coords[:, 3:] * self.translate[np.newaxis, :3])
        m[:, 3] = coords[:, 3]
        return m
433995
def cartesian_to_spherical(cartesian):
    """"""
    xyz = cartesian
    xy = xyz[:,0]**2 + xyz[:,1]**2
    r = np.sqrt(xy + xyz[:,2]**2)
    phi = np.arctan2(np.sqrt(xy), xyz[:,2]) # for elevation angle defined from Z-axis down
    #ptsnew[:,4] = np.arctan2(xyz[:,2], np.sqrt(xy)) # for elevation angle defined from XY-plane up
    theta = np.arctan2(xyz[:,1], xyz[:,0])
    return np.array([r, theta, phi]).T
1322478
def screen_cv2(self):
        """"""
        pil_image = self.screen.convert('RGB')
        cv2_image = np.array(pil_image)
        pil_image.close()
        # Convert RGB to BGR 
        cv2_image = cv2_image[:, :, ::-1]
        return cv2_image
219600
def _cartesian_to_spherical(cls, coord, center):
        
        x, y, z, vx, vy, vz = coord
        r = np.linalg.norm(coord[:3])
        phi = arcsin(z / r)
        theta = arctan2(y, x)

        r_dot = (x * vx + y * vy + z * vz) / r
        phi_dot = (vz * (x ** 2 + y ** 2) - z * (x * vx + y * vy)) / (r ** 2 * sqrt(x ** 2 + y ** 2))
        theta_dot = (x * vy - y * vx) / (x ** 2 + y ** 2)

        return np.array([r, theta, phi, r_dot, theta_dot, phi_dot], dtype=float)
694741
def ang2pix(nside, theta, phi, nest=False, lonlat=False):
    """"""
    lon, lat = _healpy_to_lonlat(theta, phi, lonlat=lonlat)
    return lonlat_to_healpix(lon, lat, nside, order='nested' if nest else 'ring')
1374539
def _agg_cv_result(raw_results, eval_train_metric=False):
    """"""
    cvmap = collections.defaultdict(list)
    metric_type = {}
    for one_result in raw_results:
        for one_line in one_result:
            if eval_train_metric:
                key = "{} {}".format(one_line[0], one_line[1])
            else:
                key = one_line[1]
            metric_type[key] = one_line[3]
            cvmap[key].append(one_line[2])
    return [('cv_agg', k, np.mean(v), metric_type[k], np.std(v)) for k, v in cvmap.items()]
1236727
def minimum_image(self):
        """"""
        if self.box_vectors is None:
            raise ValueError('No periodic vectors defined')
        else:
            self.r_array = minimum_image(self.r_array, self.box_vectors.diagonal())
        
        return self
1282249
def surface_is_valid3(surface):
    """"""
    if NO_IMAGES:
        return

    edge1, edge2, edge3 = surface.edges
    N = 128
    # Compute points on each edge.
    std_s = np.linspace(0.0, 1.0, N + 1)
    points1 = edge1.evaluate_multi(std_s)
    points2 = edge2.evaluate_multi(std_s)
    points3 = edge3.evaluate_multi(std_s)
    # Compute the actual boundary where the Jacobian is 0.
    s_vals = np.linspace(0.0, 0.2, N)
    t_discrim = np.sqrt((1.0 - s_vals) * (1.0 - 5.0 * s_vals))
    t_top = 0.5 * (1.0 - s_vals + t_discrim)
    t_bottom = 0.5 * (1.0 - s_vals - t_discrim)
    jacobian_zero_params = np.zeros((2 * N - 1, 2), order="F")
    jacobian_zero_params[:N, 0] = s_vals
    jacobian_zero_params[:N, 1] = t_top
    jacobian_zero_params[N:, 0] = s_vals[-2::-1]
    jacobian_zero_params[N:, 1] = t_bottom[-2::-1]
    jac_edge = surface.evaluate_cartesian_multi(jacobian_zero_params)
    # Add the surface to the plot and add a dashed line
    # for each "true" edge.
    figure = plt.figure()
    ax = figure.gca()
    line, = ax.plot(jac_edge[0, :], jac_edge[1, :])
    color = line.get_color()
    ax.plot(points1[0, :], points1[1, :], color="black", linestyle="dashed")
    ax.plot(points2[0, :], points2[1, :], color="black", linestyle="dashed")
    ax.plot(points3[0, :], points3[1, :], color="black", linestyle="dashed")
    polygon = np.hstack([points1[:, 1:], points2[:, 1:], jac_edge[:, 1:]])
    add_patch(ax, polygon, color, with_nodes=False)
    ax.axis("scaled")
    ax.set_xlim(-0.0625, 1.0625)
    ax.set_ylim(-0.0625, 1.0625)
    save_image(ax.figure, "surface_is_valid3.png")