2017-04-23 78 views
1

我想从经纬度边界定义的用户指定的(非矩形)域中的每个网格单元提取数据。我的输入文件位于曲线网格上。我已经尝试过各种方法python,cdo,ncks,但我仍然无法弄清楚。我只想为输入ncfile的多边形域子集中的每个网格单元输入时间序列信息。我NCFILE信息输入在这里给出:从域ncfile提取数据

$ ncdump -h 1979_sfc_out.nc 

netcdf \1979_sfc_out { 
dimensions: 
    x = 83 ; 
    y = 94 ; 
    time = UNLIMITED ; // (8736 currently) 
    nv4 = 4 ; 
variables: 
    float time(time) ; 
     time:axis = "T" ; 
     time:long_name = "time" ; 
     time:standard_name = "time" ; 
     time:units = "hours since 1979-1-2 00:00:00" ; 
     time:calendar = "standard" ; 
    float x(x) ; 
     x:axis = "x" ; 
     x:long_name = "X-coordinate in Cartesian system" ; 
     x:standard_name = "projection_x_coordinate" ; 
     x:units = "meters" ; 
    float y(y) ; 
     y:axis = "y" ; 
     y:long_name = "Y-coordinate in Cartesian system" ; 
     y:standard_name = "projection_y_coordinate" ; 
     y:units = "meters" ; 
    float lon(y, x) ; 
     lon:units = "degrees_east" ; 
     lon:valid_range = -180., 180. ; 
     lon:standard_name = "longitude" ; 
     lon:bounds = "lon_bnds" ; 
    float lat(y, x) ; 
     lat:units = "degrees_north" ; 
     lat:valid_range = -90., 90. ; 
     lat:standard_name = "latitude" ; 
     lat:bounds = "lat_bnds" ; 
    float lon_bnds(y, x, nv4) ; 
     lon_bnds:units = "degreesE" ; 
    float lat_bnds(y, x, nv4) ; 
     lat_bnds:units = "degreesN" ; 
    char mapping ; 
     mapping:false_easting = 0. ; 
     mapping:false_northing = 0. ; 
     mapping:grid_mapping_name = "polar_stereographic" ; 
     mapping:latitude_of_projection_origin = 90. ; 
     mapping:standard_parallel = 64. ; 
     mapping:straight_vertical_longitude_from_pole = -152. ; 
     mapping:semi_major_axis = 6370000. ; 
     mapping:semi_minor_axis = 6370000. ; 
    float SEAICE(time, y, x) ; 
     SEAICE:_FillValue = -9999.f ; 
     SEAICE:units = "fraction" ; 
     SEAICE:long_name = "Ice concentration (ice=1;no ice=0)" ; 
     SEAICE:grid_mapping = "mapping" ; 
     SEAICE:coordinates = "lon lat" ; 

一些我已经尝试的事情是

lat1=71.2 
lat2=72.9 
lon1=-176.5 
lon2=-160 

cdo sellonlatbox,lon1,lon2,lat1,lat2 $ifile $box1_ofile 
cdo sellonlatbox (Abort): Float parameter >lon1< contains invalid character at position 1! 

我认为这个问题是我的输入文件具有X,Y尺寸以米为单位(不有一个负号,可能是'位置1的无效字符'),并且我要求cdo以度为单位提取纬度/经度尺寸。我有我的nc文件中的变量'映射',这可能有助于将米转换为纬度/经度,但我无法弄清楚如何去做。

而且ncks在我这里不起作用,可能是由于相同的电表< - > lon/lat问题,我用cdo看到。

ncks -v SEAICE,U10,V10 -d latitude,71.2,72.9 -d longitude,-176.5,-160. $ifile -O $box1_ofile 

ncks: ERROR dimension latitude is not in input file 

虽然我想提取的多边形,这些例子我曾尝试都只是矩形子集,我认为得到多边形我可以做多个矩形子集来实现我的最终多边形形状,但如果有更好的方式做到这一点,任何意见表示赞赏。

由于

回答

1

NCO的辅助坐标feature被设计成hyperslab上非结构化网格曲线坐标与1- d纬度和经度。试试这个:

ncks -X lon_min,lon_max,lat_min,lat_max in.nc out.nc 
ncks -X -176.5,-160.,71.2,72.9 in.nc out.nc 

您也可以菊花链的多个-X选项来获得怪异的多边形。不幸的是,这可能对你不起作用,因为我刚刚注意到,你有一个2D lat和lon的曲线网格。对于那个尝试ncap2,其中feature。本手册中的示例显示了如何沿矩形边界进行遮罩,并且可以将where()语句中的条件进行菊花链连接以获取多边形。这不会更改输出文件的维度,但它允许您将面外的所有内容都设置为_FillValue。