Convert map coordinates to pixel coordinates based on geotransform
    
    Accepts float or NumPy arrays

    GDAL model used here - upper left corner of upper left pixel for mX, mY (and in GeoTransform)
973076
def mapToPixel(mX, mY, geoTransform):
    
    mX = np.asarray(mX)
    mY = np.asarray(mY)
    if geoTransform[2] + geoTransform[4] == 0:
        pX = ((mX - geoTransform[0]) / geoTransform[1]) - 0.5
        pY = ((mY - geoTransform[3]) / geoTransform[5]) - 0.5
    else:
        pX, pY = applyGeoTransform(mX, mY, invertGeoTransform(geoTransform))
    #return int(pX), int(pY)
    return pX, pY
548858
def image_corner(self, corner):
        """"""
        if corner not in self.corner_types():
            raise GeoRaster2Error('corner %s invalid, expected: %s' % (corner, self.corner_types()))

        x = 0 if corner[1] == 'l' else self.width
        y = 0 if corner[0] == 'u' else self.height
        return Point(x, y)
989196
def pixel_to_geo(pixel, level):
        """"""
        pixel_x = pixel[0]
        pixel_y = pixel[1]
        map_size = float(TileSystem.map_size(level))
        x = (TileSystem.clip(pixel_x, (0, map_size - 1)) / map_size) - 0.5
        y = 0.5 - (TileSystem.clip(pixel_y, (0, map_size - 1)) / map_size)
        lat = 90 - 360 * atan(exp(-y * 2 * pi)) / pi
        lon = 360 * x
        return round(lat, 6), round(lon, 6)
1093012
def _transform_col(self, x, i):
        
        return x.fillna(NAN_INT).map(self.label_encoders[i]).fillna(0)
334473
def get(self, x, y):
        
        x = normalize(x)
        y = normalize(y)
        dot_index = pixel_map[y % 4][x % 2]
        col, row = get_pos(x, y)
        char = self.chars.get(row, {}).get(col)

        if not char:
            return False

        if type(char) != int:
            return True

        return bool(char & dot_index)
1252277
def ds_geom(ds, t_srs=None):
    
    gt = ds.GetGeoTransform()
    ds_srs = get_ds_srs(ds)
    if t_srs is None:
        t_srs = ds_srs
    ns = ds.RasterXSize
    nl = ds.RasterYSize
    x = np.array([0, ns, ns, 0, 0], dtype=float)
    y = np.array([0, 0, nl, nl, 0], dtype=float)
    #Note: pixelToMap adds 0.5 to input coords, need to account for this here
    x -= 0.5
    y -= 0.5
    mx, my = pixelToMap(x, y, gt)
    geom_wkt = 'POLYGON(({0}))'.format(', '.join(['{0} {1}'.format(*a) for a in zip(mx,my)]))
    geom = ogr.CreateGeometryFromWkt(geom_wkt)
    geom.AssignSpatialReference(ds_srs)
    if not ds_srs.IsSame(t_srs):
        geom_transform(geom, t_srs)
    return geom
424421
def _get_y_axis(self):
        """"""
        # 0,0 is top/left top top/left pixel. Actual x/y coord of that pixel are (.5,.5).
        y_centers = np.linspace(.5, self.y_size - .5, self.y_size)
        x_centers = y_centers * 0
        return (self.geotransform[3]
                + self.geotransform[4] * x_centers
                + self.geotransform[5] * y_centers)
1078741
function (x, y, width, height) {
        var rect = this.getBoundingRect();
        var rawTransform = this._rawTransformable;

        rawTransform.transform = rect.calculateTransform(
            new BoundingRect(x, y, width, height)
        );

        rawTransform.decomposeTransform();

        this._updateTransform();
    }
1364648
def from_geo(geo, level):
    
    pixel = TileSystem.geo_to_pixel(geo, level)
    tile = TileSystem.pixel_to_tile(pixel)
    key = TileSystem.tile_to_quadkey(tile, level)
    return QuadKey(key)
405274
def get_xy_1D(ds, stride=1, getval=False):
    
    gt = ds.GetGeoTransform()
    #stride = stride_m/gt[1]
    pX = np.arange(0, ds.RasterXSize, stride)
    pY = np.arange(0, ds.RasterYSize, stride)
    mX, dummy = pixelToMap(pX, pY[0], gt)
    dummy, mY = pixelToMap(pX[0], pY, gt)
    return mX, mY
106980
def pixel(bounding_box)
      x, y = MapKit.latlng2pixel(@lat, @lng, bounding_box.zoom)
      tile_x, tile_y = MapKit.latlng2pixel(bounding_box.top, bounding_box.left, bounding_box.zoom)
      [x-tile_x, y-tile_y]
    end
604532
def index_to_coordinates(self, index):
        
        self.validate_index(index)
        x = int(index % self.length)
        y = int((index - x) / self.length)
        return x, y
905669
def reproject_to_grid_coordinates(self, grid_coordinates, interp=gdalconst.GRA_NearestNeighbour):
        
        source_dataset = self.grid_coordinates._as_gdal_dataset()
        dest_dataset = grid_coordinates._as_gdal_dataset()
        rb = source_dataset.GetRasterBand(1)
        rb.SetNoDataValue(NO_DATA_VALUE)
        rb.WriteArray(np.ma.filled(self.raster_data, NO_DATA_VALUE))

        gdal.ReprojectImage(source_dataset, dest_dataset,
                            source_dataset.GetProjection(),
                            dest_dataset.GetProjection(),
                            interp)
        dest_layer = self.clone_traits()
        dest_layer.grid_coordinates = grid_coordinates
        rb = dest_dataset.GetRasterBand(1)
        dest_layer.raster_data = np.ma.masked_values(rb.ReadAsArray(), NO_DATA_VALUE)
        return dest_layer
1310374
def ds_cT(ds, x, y, xy_srs=wgs_srs):
    
    #Convert lat/lon to projected srs
    ds_srs = get_ds_srs(ds)
    #If xy_srs is undefined, assume it is the same as ds_srs
    mX = x
    mY = y
    if xy_srs is not None:
        if not ds_srs.IsSame(xy_srs):
            mX, mY, mZ = cT_helper(x, y, 0, xy_srs, ds_srs)
    return mX, mY
1119480
def _load_GeoTransform(self):
        """"""
        def load_lon():
            return arange(ds.RasterXSize)*b[1]+b[0]

        def load_lat():
            return arange(ds.RasterYSize)*b[5]+b[3]
        ds = self.ds
        b = self.ds.GetGeoTransform()  # bbox, interval
        if with_dask:
            lat = Array(
                {('lat', 0): (load_lat,)}, 'lat', (self.ds.RasterYSize,),
                shape=(self.ds.RasterYSize,), dtype=float)
            lon = Array(
                {('lon', 0): (load_lon,)}, 'lon', (self.ds.RasterXSize,),
                shape=(self.ds.RasterXSize,), dtype=float)
        else:
            lat = load_lat()
            lon = load_lon()
        return Variable(('lat',), lat), Variable(('lon',), lon)
967455
def warp_image_by_corner_points_projection(corner_points, image):
    
    # Clarify by storing in named variables.
    top_left, top_right, bottom_left, bottom_right = np.array(corner_points)

    top_edge = np.linalg.norm(top_right - top_left)
    bottom_edge = np.linalg.norm(bottom_right - bottom_left)
    left_edge = np.linalg.norm(top_left - bottom_left)
    right_edge = np.linalg.norm(top_right - bottom_right)

    L = int(np.ceil(max([top_edge, bottom_edge, left_edge, right_edge])))
    src = np.array([top_left, top_right, bottom_left, bottom_right])
    dst = np.array([[0, 0], [L - 1, 0], [0, L - 1], [L - 1, L - 1]])

    tr = ProjectiveTransform()
    tr.estimate(dst, src)
    warped_image = warp(image, tr, output_shape=(L, L))
    out = resize(warped_image, (500, 500))

    return out
576197
def from_scanner(self, x, y, z):
        
        if self.transform is None:
            raise ValueError("No transform set for MRSData object {}".format(self))

        transformed_point = numpy.linalg.inv(self.transform) * numpy.matrix([x, y, z, 1]).T

        return numpy.squeeze(numpy.asarray(transformed_point))[0:3]