如何将Lambert Conic Conformal栅格投影转换为经纬度度数R

Question

如何将Lambert Conic Conformal栅格投影转换为经纬度度数R

3

我有一个栅格图像，从netcdf文件中获取，它采用兰伯特等角投影。

    library(meteoForecast)
    wrf_temporary <- getRaster("temp", day = Sys.Date(), frames = 'complete', resolution = 36, service = "meteogalicia")
    wrf_temporary

extent      : -18, 4230, -18, 3726  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=lcc +lat_1=43 +lat_2=43 +lat_0=34.82300186157227 +lon_0=-14.10000038146973 +x_0=536402.34 +y_0=-18558.61 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=km +no_defs

现在我想把那个wrf_temporary栅格转换为"+proj=longlat +datum=WGS84"(纬度经度度数)。应该怎么做呢？我希望得到这样的结果：

    mfExtent('meteogalicia', resolution = 36)

class       : Extent 
xmin        : -49.18259 
xmax        : 18.789 
ymin        : 24.03791 
ymax        : 56.06608

已经尝试了很多选项，但没有一个能给出正确的结果...

- Ricardo Faria

1

你尝试过使用projectRaster吗？如果有，你是怎么做的？如果没有，可以试试看。 - Spacedman

是的，尝试过 projectRaster(wrf_temporary, crs="+proj=longlat +datum=WGS84")，但没有成功 :( - Ricardo Faria

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Spacedman · Accepted Answer

这应该可以工作：

> ll = projectRaster(wrf_temporary,crs="+init=epsg:4326")
> plot(ll[[1]])

并生成了这个:

乍一看似乎没问题，但是格林威治子午线（经度=0）并没有穿过伦敦。

当分辨率设置为其他值时，4或12，就不会出现这个问题。当以resolution=36调用时，你会得到一个比其他raster更大的范围，但具有相同的投影字符串。这肯定是不对的。例如，以下图中的（0,0）坐标应该是地球上的同一点，因为它们声称是相同的投影，但实际上它们不是。

我认为分辨率为12的那些是正确的。下面是将该raster转换为经纬度，并覆盖EU矢量覆盖层的结果：

完全匹配。

所以我认为这是一个bug，要么getRaster猜测投影并且出错了，要么在裁剪后没有应用投影变换，或者使用的任何服务都在谎称投影。

getRaster正在猜测。下载的NetCDF文件中包含投影信息，但是以另一种格式呈现。对于resolution=36文件，投影信息部分如下：

int Lambert_Conformal ;
    Lambert_Conformal:standard_parallel = 43., 43. ;
    Lambert_Conformal:latitude_of_projection_origin = 24.2280006408691 ;
    Lambert_Conformal:false_easting = 2182.62935 ;
    Lambert_Conformal:false_northing = -269.65597 ;
    Lambert_Conformal:grid_mapping_name = "lambert_conformal_conic" ;
    Lambert_Conformal:longitude_of_central_meridian = -14.1000003814697 ;
    Lambert_Conformal:_CoordinateTransformType = "Projection" ;
    Lambert_Conformal:_CoordinateAxisTypes = "GeoX GeoY" ;

这与R包使用的resolution=12投影有些不同。因此，请制作一个符合要求的投影：

p = "+proj=lcc +lat_1=43 +lat_2=43 +lat_0=24.2280006408691 +lon_0=-14.1000003814697 +x_0=2182629.35 +y_0=-269655.97 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=km +no_defs"
 wrf <- getRaster('temp', day = testDay, resolution=36)
 projection(wrf) = p

那么，带有欧盟叠加层的测试...

请注意，这次我重新投影了欧盟。对经纬度进行projectRaster也可以：

> wrfLL = projectRaster(wrf, crs="+init=epsg:4326")
> plot(wrfLL[[1]])
> abline(v=0)

伦敦子午线恢复原位。

注：伦敦子午线是指通过英国伦敦郊区的格林尼治天文台的经线，被认为是世界标准时间的起点。