NCL 对nc文件进行区域掩码

haysen

就是比如有全国的降水数据，想从中获取其中某个区域的降水数据（比如只获取黄河流域的），其他的全部制成无效值

需要用到shapefile_utils.ncl  这个文件中的函数，这个文件在官网上能下到，需要流域shape文件

1. load "$NCARG_ROOT/shapefile_utils.ncl"
   
2. 
   
3. begin
   
4. 
   
5. f      = addfile("prcp.nc", "r")      ;//打开NC文件
   
6. prcp_temp    = f->prcp(0,:,:)                      ;//读取数据
   
7. 
   
8. shpfile   = "./liuyu.shp"
   
9. print_shapefile_info(shpfile)
   
10. plot_shapefile(shpfile)
    
11. 
    
12. 
    
13. opt             = True
    
14. opt@return_mask = True
    
15. opt@shape_var   = "W1102WB0_1"
    
16. opt@shape_names = "黄河流域片"
    
17. 
    
18. land_mask  = shapefile_mask_data(prcp_temp,shpfile,opt)
    
19. 
    
20. prcp_land_mask  =  where(land_mask.eq.1,prcp_temp,prcp_temp@_FillValue)
    
21. 
    
22. copy_VarMeta(prcp_temp,prcp_land_mask)
    
23. 
    
24. pcp = prcp_land_mask
    
25. 
    
26. end

复制代码

yysunnyboy

谢谢大神的经验分享！   

qingshan

楼主好，我按您的代码运行出现如下错误，求指教
BDQS@DESKTOP-97QP3T6 ~
$ ncl QTPmask.ncl
Copyright (C) 1995-2019 - All Rights Reserved
University Corporation for Atmospheric Research
NCAR Command Language Version 6.6.2
The use of this software is governed by a License Agreement.
See http://www.ncl.ucar.edu/ for more details.
warning:_NclOpenFile: Can not open file <./DBATP_Polygon.shp>; file format not supported or file is corrupted
warning:getfilevarnames: f is not a valid file variable
fatal:file (f) isn't defined
fatal:["Execute.c":8637]:Execute: Error occurred at or near line 48 in file ./shapefile_utils.ncl
fatal:["Execute.c":8637]:Execute: Error occurred at or near line 11 in file QTPmask.ncl

以下是我的代码
;load "./gsn_code.ncl"
;load "./gsn_csm.ncl"
;load "./contributed.ncl"
load "./shapefile_utils.ncl"
begin

f      = addfile("1979.nc", "r")      ;//打开NC文件
prcp_temp    = f->temp(0,:,:)                      ;//读取数据

shpfile   = "./DBATP_Polygon.shp"
print_shapefile_info(shpfile)
plot_shapefile(shpfile)


opt             = True
opt@return_mask = True
opt@shape_var   = "W1102WB0_1"
opt@shape_names = "QTP"

land_mask  = shapefile_mask_data(prcp_temp,shpfile,opt)

prcp_land_mask  =  where(land_mask.eq.1,prcp_temp,prcp_temp@_FillValue)

copy_VarMeta(prcp_temp,prcp_land_mask)

pcp = prcp_land_mask

end  

会心一笑

qingshan 发表于 2020-3-13 13:30
楼主好，我按您的代码运行出现如下错误，求指教
BDQS@DESKTOP-97QP3T6 ~
$ ncl QTPmask.ncl

我也是出现这种错误，同问，谢谢