NLLGrid.py 文件源码

python
阅读 20 收藏 0 点赞 0 评论 0

项目:nllgrid 作者: claudiodsf 项目源码 文件源码
def project(self, lon, lat):
        """Project lon, lat into grid coordinates."""
        if self.proj_name == 'LAMBERT':
            ellipsoids = {
                'WGS-84': 'WGS84',
                'GRS-80': 'GRS80',
                'WGS-72': 'WGS72',
                'Australian': 'aust_SA',
                'Krasovsky': 'krass',
                'International': 'new_intl',
                'Hayford-1909': 'intl',
                'Clarke-1880': 'clrk80',
                'Clarke-1866': 'clrk66',
                'Airy': 'airy',
                'Bessel': 'bessel',
                # 'Hayford-1830':
                'Sphere': 'sphere'
                }
            try:
                ellps = ellipsoids[self.proj_ellipsoid]
            except KeyError:
                raise ValueError(
                    'Ellipsoid not supported: {}'.format(self.proj_ellipsoid))
            p = Proj(proj='lcc', lat_0=self.orig_lat, lon_0=self.orig_lon,
                     lat_1=self.first_std_paral, lat_2=self.second_std_paral,
                     ellps=ellps)
        elif self.proj_name == 'SIMPLE':
            p = Proj(proj='eqc', lat_0=self.orig_lat, lon_0=self.orig_lon)
        else:
            raise ValueError('Projection not supported: {}'.format(
                self.proj_name))
        x, y = p(lon, lat)
        x = np.array(x)
        y = np.array(y)
        x[np.isnan(lon)] = np.nan
        y[np.isnan(lat)] = np.nan
        x /= 1000.
        y /= 1000.
        # inverse rotation
        theta = np.radians(self.map_rot)
        x1 = x*np.cos(theta) + y*np.sin(theta)
        y1 = -x*np.sin(theta) + y*np.cos(theta)
        x = x1
        y = y1
        # Try to return the same type of lon, lat
        if not isinstance(lon, np.ndarray):
            if isinstance(lon, Iterable):
                x = type(lon)(x)
            else:
                x = float(x)
        if not isinstance(lat, np.ndarray):
            if isinstance(lat, Iterable):
                y = type(lat)(y)
            else:
                y = float(y)
        return x, y
评论列表
文章目录


问题


面经


文章

微信
公众号

扫码关注公众号