Failure to create oceanfile

Spatial Allocator
1

Hi~

I am using spatial allocator to create ocean file for CMAQ. The domain is out of North America, so I made the shape files myself and import them to Spatial Allocator. However Spatial Allocator stopped without any error information. Here is my log file

WARNING: Environment variable: MAX_INPUT_FILE_SHAPES, not set
WARNING: Environment variable: OVERLAY_TYPE, not set
WARNING: Environment variable: MAX_LINE_SEG, not set
WARNING: Environment variable: MAX_LINE_SEG, not set
WARNING: polyIsect: There is no intersection between the sets of polygons
WARNING: Possible empty intersection in data polygons
Spatial Allocator Version 3.6, 03/10/2009

EV: MIMS_PROCESSING=ALLOCATE
EV: ALLOCATE_ATTRS=TYPE
EV: OUTPUT_FILE_NAME=/home/xiao/Build_WRF/Spatial-Allocator/output/ocean_file_Grid-days20.ncf
Using ALLOCATE mode
EV: OUTPUT_FILE_TYPE=IoapiFile
EV: OUTPUT_FILE_TYPE=IoapiFile
EV: OUTPUT_FILE_TYPE=IoapiFile
Reading Regular Grid

EV: OUTPUT_GRID_NAME=Grid-days20
MAX_LINE_SEG not set, discretization intervals disabled
griddesc file name = /home/xiao/Build_WRF/Spatial-Allocator/data/GRIDDESC

Ellipsoid var = OUTPUT_FILE_ELLIPSOID
EV: OUTPUT_FILE_ELLIPSOID=+a=6370997.0,+b=6370997.0
Ellipsoid=+a=6370997.0,+b=6370997.0
EV: OUTPUT_GRID_NAME=Grid-days20
Not using BB optimization

Finished reading the output fileā€¦

EV: INPUT_FILE_TYPE=ShapeFile
EV: INPUT_FILE_ELLIPSOID=+a=6370997.0,+b=6370997.0
Ellipsoid=+a=6370997.0,+b=6370997.0
EV: INPUT_FILE_MAP_PRJN=+proj=lcc,+lat_1=30,+lat_2=60,+lat_0=40,+lon_0=112
EV: INPUT_FILE_TYPE=ShapeFile
EV: INPUT_FILE_NAME=/home/xiao/Build_WRF/Spatial-Allocator/data/surfzone/ChinaMapCoast
MAX_LINE_SEG not set, discretization intervals disabled
max_line_seg=

Reading Shapefile /home/xiao/Build_WRF/Spatial-Allocator/data/surfzone/ChinaMapCoast
Shapefile Type: 5 Polygon # of Shapes: 44

Input projection:

args[1]=+a=6370997.0

args[2]=+b=6370997.0

param 0 = +proj=lcc,
param 3 = +lat_1=30,
param 4 = +lat_2=60,
param 5 = +lat_0=40,
param 6 = +lon_0=112,
PROJ args=+proj=lcc+a=6370997.0+b=6370997.0+lat_1=30+lat_2=60+lat_0=40+lon_0=112

Output projection:

args[1]=+a=6370997.0

args[2]=+b=6370997.0

setting Lambert Conic Conformal ā€“ 2SP parameters
PROJ args=+proj=lcc+a=6370997.0+b=6370997.0+lat_1=30.000000+lat_2=60.000000+lon_0=112.000000+lat_0=40.000000+units=m

Total area = 9.9707e+12

Skipped 0 polygons and 0 vertices

File bounding box:
xmin = 16204787.444 xmax = 21277082.714 ymin = 2009033.435 ymax = 6010188.757

Finished reading the input shape fileā€¦

EV: ALLOCATE_ATTRS=TYPE
EV: ALLOCATE_ATTRS=TYPE
Not using function for weights

Using attribute TYPE

Finished reading the input attribute fileā€¦

xmin = 16204787.444 xmax = 21277082.714 ymin = 2009033.435 ymax = 6010188.757

xmin = -789000.000 xmax = 1367000.000 ymin = -2165500.000 ymax = -9500.000

0.0u 0.0s 0:00.15 46.6% 0+0k 0+8io 0pf+0w

Does anyone know the reason? Thank you in advance!

2

If the spatial allocator script is set up to write the OUTPUT to a directory, and that directory doesnā€™t exist, then the spatial allocator will fail to create the output file. Please add a check to the script that if the directory doesnā€™t exist to make it prior to the run.

setenv OUTDIR $cwd/output

if ( ! -d ā€œ$OUTDIRā€ ) mkdir -p $OUTDIR

3

Hi Liz,

Thanks for your reply! I wrote the code to the script. However the log file does not change. Here is my ā€˜alloc_surf_zone_to_oceanfile_China.cshā€™ script

cd /home/xiao/Build_WRF/Spatial-Allocator/bin
source sa_setup.csh
cd ā€¦/scripts

setenv DEBUG_OUTPUT Y

Set executable

setenv EXE ā€œ$SA_HOME/bin/64bits/allocator.exeā€

Set Input Directory

setenv DATADIR $SA_HOME/data
setenv OUTPUT $SA_HOME/output

if ( ! -d ā€œ$OUTPUTā€ ) mkdir -p $OUTPUT

Select method of spatial analysis

setenv MIMS_PROCESSING ALLOCATE

setenv TIME time

#set ā€œdataā€ shapefile parameters
setenv GRIDDESC $DATADIR/GRIDDESC

#set parameters for file being allocated
setenv INPUT_FILE_NAME $DATADIR/surfzone/ChinaMapCoast
setenv INPUT_FILE_TYPE ShapeFile
setenv INPUT_FILE_MAP_PRJN ā€œ+proj=lcc,+lat_1=30,+lat_2=60,+lat_0=40,+lon_0=112ā€
setenv INPUT_FILE_ELLIPSOID ā€œ+a=6370997.0,+b=6370997.0ā€
setenv ALLOCATE_ATTRS TYPE
setenv ALLOC_MODE_FILE ALL_AREAPERCENT

#Set this to SURF_ZONE to create the variables needed for CMAQ OCEANfile
setenv ALLOC_ATTR_TYPE SURF_ZONE

Set name and path of resulting shapefile

setenv OUTPUT_FILE_TYPE IoapiFile
setenv OUTPUT_GRID_NAME Grid-days20
setenv OUTPUT_FILE_MAP_PRJN ā€œ+proj=lcc,+lat_1=30,+lat_2=60,+lat_0=40,+lon_0=112ā€
setenv OUTPUT_FILE_ELLIPSOID ā€œ+a=6370997.0,+b=6370997.0ā€
setenv OUTPUT_FILE_NAME $OUTPUT/ocean_file.ncf

#echo ā€œAllocating surf zone data to CMAQ OCEANfileā€
$TIME $EXE

I guess I have the output directory. Could you please take a look at it? Thanks!

4

From the warning message, it means that your ocean surfzone shape file doesnā€™t intersect with your grid domain. Please double check if the projection your shape file used is same as what you specified in script. Also check your GRIDDESC on the projection used for Grid-days20. If possible, your can do some QA using ArcGIS.

1 Like
5

Hi Yang,

Thanks! And you are right. The projections of the shapefile and GRIDDESC are not aligned. I correct the projection of the shapefile with QGIS. Now everything goes well.