Python pyproj.Proj方法代码示例

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def dbcoordinates(db, pier, epsgcode='epsg:4326'):

    try:
        from pyproj import Proj, transform
    except ImportError:
        print ("dbcoordinates: You need to install pyproj to use this method")
        return (0.0 , 0.0)

    startlong = dbselect(db,'PierLong','PIERS','PierID = "A2"')
    startlat = dbselect(db,'PierLat','PIERS','PierID = "A2"')
    coordsys = dbselect(db,'PierCoordinateSystem','PIERS','PierID = "A2"')
    startlong = float(startlong[0].replace(',','.'))
    startlat = float(startlat[0].replace(',','.'))
    coordsys = coordsys[0].split(',')[1].lower().replace(' ','')

    # projection 1: GK M34
    p1 = Proj(init=coordsys)
    # projection 2: WGS 84
    p2 = Proj(init='epsg:4326')
    # transform this point to projection 2 coordinates.
    lon1, lat1 = transform(p1,p2,startlong,startlat)

    return (lon1,lat1) 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def GetUTMEastingNorthing(EPSG_string,latitude,longitude):
    """This returns the easting and northing for a given latitude and longitide

    Args:
        ESPG_string (str): The ESPG code. 326XX is for UTM north and 327XX is for UTM south
        latitude (float): The latitude in WGS84
        longitude (float): The longitude in WGS84

    Returns:
        easting,northing The easting and northing in the UTM zone of your selection

    Author:
        Simon M Mudd
    """

    #print "Yo, getting this stuff: "+EPSG_string
    # The lat long are in epsg 4326 which is WGS84
    inProj = Proj(init='epsg:4326')
    outProj = Proj(init=EPSG_string)
    ea,no = transform(inProj,outProj,longitude,latitude)

    return ea,no 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def get_unit_multiplier_from_epsg(epsg):
    # Can't find any lightweight way to do this, but I depend heavily on pyproj right now
    # Until there's a better way, get the unit by converting (1,0) from 'unit' to 'meter'
    meter_proj = pyproj.Proj(init='epsg:'+str(epsg), preserve_units=False) # Forces meter
    unit_proj = pyproj.Proj(init='epsg:'+str(epsg), preserve_units=True) # Stays in native unit
    try:
        # Due to numerical precision and the requirements for foot vs survey foot
        # Use a larger number for the transform
        scale_value = 1.0e6
        x2, y2 = pyproj.transform(unit_proj, meter_proj, scale_value, 0.0)
        return x2/scale_value
    except:
        pass
    return 0.0

# Some of the epsgs found in lidar files don't represent projected coordinate systems
# but are instead the datum or other geographic reference systems 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def __init__(self, shape=None, mask=None, nodata=None,
                 crs=pyproj.Proj(_pyproj_init), y_coord_ix=0, x_coord_ix=1):
        if shape is not None:
            self.shape = shape
        else:
            self.shape = (0,0)
        self.crs = crs
        if nodata is None:
            self.nodata = np.nan
        else:
            self.nodata = nodata
        if mask is None:
            self.mask = np.ones(shape).astype(bool)
        else:
            self.mask = mask
        self.y_coord_ix = y_coord_ix
        self.x_coord_ix = x_coord_ix 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def test_projplot(self):

        # this caused many problems on Linux Mint distributions
        # this is just to be sure that on your system, everything is fine
        import pyproj
        import matplotlib.pyplot as plt

        pyproj.Proj(proj='latlong', datum='WGS84')
        plt.figure()
        plt.close()

        srs = ('+units=m +proj=lcc +lat_1=29.0 +lat_2=29.0 '
               '+lat_0=29.0 +lon_0=89.8')

        proj_out = pyproj.Proj("EPSG:4326", preserve_units=True)
        proj_in = pyproj.Proj(srs, preserve_units=True)

        from salem.gis import transform_proj
        lon, lat = transform_proj(proj_in, proj_out, -2235000, -2235000)
        np.testing.assert_allclose(lon, 70.75731, atol=1e-5) 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def transform(p1, p2, c1, c2, c3):
    if PYPROJ:
        if type(p1) is Proj and type(p2) is Proj:
            if p1.srs != p2.srs:
                return proj_transform(p1, p2, c1, c2, c3)
            else:
                return (c1, c2, c3)
        elif type(p2) is TransverseMercator:
            wgs84 = Proj(init="EPSG:4326")
            if p1.srs != wgs84.srs:
                t2, t1, t3 = proj_transform(p1, wgs84, c1, c2, c3)
            else:
                t1, t2, t3 = c2, c1, c3  # mind c2, c1 inversion
            tm1, tm2 = p2.fromGeographic(t1, t2)
            return (tm1, tm2, t3)
    else:
        if p1.spherical:
            t1, t2 = p2.fromGeographic(c2, c1)  # mind c2, c1 inversion
            return (t1, t2, c3)
        else:
            return (c1, c2, c3) 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def _transform_proj(x1, y1, epsg_src, epsg_dst, polar=False):
    if PYPROJ_VER == 1:
        proj_src = Proj(init='epsg:{0}'.format(epsg_src))
        _log_proj_crs(proj_src, proj_desc='src', espg=epsg_src)
        proj_dst = Proj(init='epsg:{0}'.format(epsg_dst))
        _log_proj_crs(proj_dst, proj_desc='dst', espg=epsg_dst)
        x2, y2 = transform(proj_src, proj_dst, x1, y1)
    elif PYPROJ_VER == 2:
        # With PROJ 6+ input axis ordering needs honored per projection, even
        #   though always_xy should fix it (doesn't seem to work for UPS)
        crs_src = CRS.from_epsg(epsg_src)
        _log_proj_crs(crs_src, proj_desc='src', espg=epsg_src)
        crs_dst = CRS.from_epsg(epsg_dst)
        _log_proj_crs(crs_dst, proj_desc='dst', espg=epsg_dst)
        ct = Transformer.from_crs(crs_src, crs_dst, always_xy=(not polar))
        if polar:
            y2, x2 = ct.transform(y1, x1)
        else:
            x2, y2 = ct.transform(x1, y1)
    else:
        raise MgrsException('pyproj version unsupported')

    return x2, y2 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def _proj(**kwargs):
    projection = kwargs.get("projection")
    if projection == "WGS84":
        proj = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
    if projection == "NAD27":
        proj = "+proj=longlat +ellps=clrk66 +datum=NAD27 +no_defs"
    if projection == "UTM":
        zone = _utm_zone(kwargs.get("longitude"))
        proj = "+proj=utm +zone={0:d} +datum=WGS84 +units=m +no_defs"
        proj = proj.format(zone)
    if projection == "LCC":
        lon0 = kwargs.get("lon0")
        lat0 = kwargs.get("lat0")
        parallel_1 = kwargs.get("parallel_1")
        parallel_2 = kwargs.get("parallel_2")
        proj = "+proj=lcc +lon_0={} +lat_0={} +lat_1={} +lat_2={} +datum=WGS84 +units=m +no_defs"
        proj = proj.format(float(lon0), float(lat0),
                           float(parallel_1), float(parallel_2))
    if projection == "TM":
        lon = kwargs.get("lon")
        lat = kwargs.get("lat")
        proj = "+proj=tmerc +lon_0={} +lat_0={} +datum=WGS84 +units=m +no_defs"
        proj = proj.format(float(lon), float(lat))

    return pyproj.Proj(proj) 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def get_area_acres(geometry):
    """ Calculate area in acres for a GeoJSON geometry
    :param geometry: A GeoJSON Polygon geometry
    :returns: Area in acres
    """

    shapely_geometry = shape(geometry)
    geom_aea = transform(
        partial(
            pyproj.transform,
            pyproj.Proj(init="EPSG:4326"),
            pyproj.Proj(
                proj="aea",
                lat1=shapely_geometry.bounds[1],
                lat2=shapely_geometry.bounds[3],
            ),
        ),
        shapely_geometry,
    )
    return round(geom_aea.area / 4046.8564224) 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def __geodesic_point_buffer(self, lon, lat, km):
        """
        Based on: https://gis.stackexchange.com/questions/289044/creating-buffer-circle-x-kilometers-from-point-using-python

        Args:
            lon:
            lat:
            km:

        Returns:
            a list of coordinates with radius km and center (lat,long) in this projection
        """
        proj_wgs84 = pyproj.Proj(init='epsg:4326')
        # Azimuthal equidistant projection
        aeqd_proj = '+proj=aeqd +lat_0={lat} +lon_0={lon} +x_0=0 +y_0=0'
        project = partial(
            pyproj.transform,
            pyproj.Proj(aeqd_proj.format(lat=lat, lon=lon)),
            proj_wgs84)
        buf = Point(0, 0).buffer(km * 1000)  # distance in metres
        return transform(project, buf).exterior.coords[:] 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def __init__(self, crs=None, epsg=None):
        if crs is not None and epsg is not None:
            raise ValueError('crs and epsg cannot both be specified')

        if epsg is not None:
            crs = _crs_from_epsg(epsg)
        if crs is not None:
            import pyproj
            crs_out = _crs_from_epsg(4326)
            proj_in = pyproj.Proj(preserve_units=True, **crs)
            proj_out = pyproj.Proj(preserve_units=True, **crs_out)
            self.transformfunc = partial(pyproj.transform, proj_in, proj_out)
        else:
            self.transformfunc = None

        self._features = [] 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def __init__(
        self,
        traffic: "Traffic",
        clustering: ClusteringProtocol,
        nb_samples: Optional[int],
        features: List[str],
        *args,
        projection: Union[None, crs.Projection, pyproj.Proj] = None,
        transform: Optional[TransformerProtocol] = None,
    ) -> None:

        self.traffic = traffic
        self.clustering = clustering
        self.nb_samples = nb_samples
        self.features = features
        self.projection = projection
        self.transform = transform
        self.X = None 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def getproj(self, withgrid=False, projformat='pyproj'):
        from PseudoNetCDF.coordutil import getproj4_from_cf_var
        gridmapping = self.variables['crs']
        proj4str = getproj4_from_cf_var(gridmapping, withgrid=withgrid)
        preserve_units = withgrid
        if projformat == 'proj4':
            return proj4str
        elif projformat == 'pyproj':
            import pyproj
            # pyproj adds +units=m, which is not right for latlon/lonlat
            if '+proj=lonlat' in proj4str or '+proj=latlon' in proj4str:
                preserve_units = True
            return pyproj.Proj(proj4str, preserve_units=preserve_units)
        elif projformat == 'wkt':
            import osr
            srs = osr.SpatialReference()
            # Imports WKT to Spatial Reference Object
            srs.ImportFromProj4(proj4str)
            srs.ExportToWkt()  # converts the WKT to an ESRI-compatible format
            return srs.ExportToWkt()
        else:
            raise ValueError('projformat must be pyproj, proj4 or wkt') 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def geodesic_point_buffer(lat, lon, meters, envelope=False):

    # get WGS 84 proj
    proj_wgs84 = pyproj.Proj(init='epsg:4326')

    # Azimuthal equidistant projection
    aeqd_proj = '+proj=aeqd +lat_0={lat} +lon_0={lon} +x_0=0 +y_0=0'
    project = partial(
        pyproj.transform,
        pyproj.Proj(aeqd_proj.format(lat=lat, lon=lon)),
        proj_wgs84)

    buf = Point(0, 0).buffer(meters)  # distance in metres

    if envelope is True:
        geom = Polygon(transform(project, buf).exterior.coords[:]).envelope
    else:
        geom = Polygon(transform(project, buf).exterior.coords[:])

    return geom.to_wkt() 

# 需要导入模块: import pyproj [as 别名]
# 或者: from pyproj import Proj [as 别名]
def projectShapes(features, toCRS):
    import pyproj
    from functools import partial
    import fiona.crs as fcrs
    from shapely.geometry import shape, mapping
    from shapely.ops import transform as shpTrans

    project = partial(
        pyproj.transform,
        pyproj.Proj(fcrs.from_epsg(4326)),
        pyproj.Proj(toCRS))

    return list(
        {'geometry': mapping(
            shpTrans(
                project,
                shape(feat['geometry']))
        )} for feat in features
        )