Project

General

Profile

Converting .shp to .nc

Added by Josh Butterworth over 5 years ago

Hi guys,

I've been reading through an old post by Karin showing how you can use NCL to convert a .shp file to .nc (below). This is to create a netcdf file of marine protected areas in the North Sea for which I have a shapefile.

;-- copy shapefile_utils.ncl from 
;-- https://www.ncl.ucar.edu/Applications/Scripts/shapefile_utils.ncl
load "$HOME/NCL/shapefiles/shapefile_utils.ncl" 

shpname  = "./India/india.shp" 
maskname = "./India/india_mask.nc" 

print_shapefile_info(shpname)

;-- open data file to get the grid to be used
f = addfile("monmean2007_celcius.nc","r")

;-- read variable
var                    =  f->T2(0,:,:)
var@lat2d              =  f->XLAT
var@lon2d              =  f->XLONG

;-- shapefile mask resources
opt             =  True
opt@return_mask =  True    ;-- this forces the return of a 0s and 1s mask array

;-- create the mask based on the given shapefile
mask_array             =  shapefile_mask_data(var, shpname, opt)
mask_array!0           = "y" 
mask_array!1           = "x" 
mask_array@coordinates = "XLAT XLONG" 

;-- create new netCDF file and write mask array
system("rm -f " + maskname)
fout = addfile(maskname,"c")

fout->mask_array =  mask_array
fout->XLAT       =  f->XLAT
fout->XLONG      =  f->XLONG

Using my data (attached) I've managed to get to the line below by substituting the T2 with time, XLAT with lat and XLONG with lon.

mask_array = shapefile_mask_data(var,shpname,opt)

But then I get the following error message:

(0)     shapefile_mask_data: Error: not a valid rectilinear, curvilinear, or uns
tructured grid

Does anyone know why this is happening and where I've gone wrong?

Many thanks,

Josh

SST.nc (204 KB) SST.nc Data with desired grid
shapefiles.tar (7.83 MB) shapefiles.tar Shapefiles

Replies (3)

RE: Converting .shp to .nc - Added by Josh Butterworth over 5 years ago

I realise in the file I uploaded after removing the time axis to enable the upload I no longer have a 'time' var, can this be done without including 'time'? If not I've attached the same data but with a years worth of time values.

SST1.nc (2.34 MB) SST1.nc grid with time axis

RE: Converting .shp to .nc - Added by Karin Meier-Fleischer over 5 years ago

Hi Josh,

this is an NCL question but I've seen that you have send it to ncl-talk too. S, you will get two answers. ;)

It doesn't matter if you have a time dimension. The adaption of the NCL script above gives a netCDF file of the masked variable.

load "$HOME/NCL/shapefiles/shapefile_utils.ncl" 

shpname  = "ospar_polygon_wdpa_simplified.shp" 
maskname = "ospar_polygon_wdpa_simplified.nc" 

print_shapefile_info(shpname)

;-- open data file to get the grid to be used
f = addfile("SST.nc","r")

;-- read variable
var                    =  f->SST(0,:,:)
var@lat1d              =  f->LAT
var@lon1d              =  f->LON

;-- shapefile mask resources
opt             =  True
opt@return_mask =  True    ;-- this forces the return of a 0s and 1s mask array

;-- create the mask based on the given shapefile
mask_array             =  shapefile_mask_data(var, shpname, opt)
mask_array!0           = "y" 
mask_array!1           = "x" 
mask_array@coordinates = "LAT LON" 

;-- create new netCDF file and write mask array
system("rm -f " + maskname)
fout = addfile(maskname,"c")

fout->mask_array =  mask_array
fout->LAT        =  f->LAT
fout->LON        =  f->LON

ncdump -h ospar_polygon_wdpa_simplified.nc

netcdf ospar_polygon_wdpa_simplified {
dimensions:
    y = 121 ;
    x = 211 ;
variables:
    int mask_array(y, x) ;
        mask_array:coordinates = "LAT LON" ;
        mask_array:_FillValue = -2147483647 ;
}
(ncl6.6) [~/Downloads] > ncdump -h ospar_polygon_wdpa_simplified.nc
netcdf ospar_polygon_wdpa_simplified {
dimensions:
    y = 121 ;
    x = 211 ;
    LAT = 121 ;
    LON = 211 ;
variables:
    int mask_array(y, x) ;
        mask_array:coordinates = "LAT LON" ;
        mask_array:_FillValue = -2147483647 ;
    float LAT(LAT) ;
        LAT:long_name = "Latitude" ;
        LAT:units = "degrees_north" ;
        LAT:axis = "Y" ;
        LAT:bounds = "LAT_bnds" ;
        LAT:point_spacing = "uneven" ;
        LAT:standard_name = "latitude" ;
    float LON(LON) ;
        LON:long_name = "Longitude" ;
        LON:units = "degrees_east" ;
        LON:axis = "X" ;
        LON:bounds = "LON_bnds" ;
        LON:modulo = 360.f ;
        LON:point_spacing = "uneven" ;
        LON:standard_name = "longitude" ;
}

-Karin

RE: Converting .shp to .nc - Added by Dana kh over 2 years ago

why the mask.nc after converting the type of lat and lon 2D not Geo2d
how to fix it ???

    (1-3/3)