[gradsusr] Generating png Image and georeferencing with GrADS
Stephen Woodbridge
stephenwoodbridge37 at gmail.com
Wed Dec 30 10:32:20 EST 2020
Followup, adding a little more information in the hope that someone will
have an idea on why my images are getting squashed.
So code for the ctl file below goes something like:
lat1=30.0
lat2=50.0
lon1=230.0
lon2=245.0
xpix=3000
ypix=4000
'reinit'
'open 'ctlfile
'set mpdset hires'
'set grads off'
'set grid off'
'set rgb 69 0 49 156'
'set background 3'
'c'
'set lat 'lat1' 'lat2
'set lon 'lon1' 'lon2
* lots of 'set rgb ... '
'set gxout barb'
'set gxout shaded'
'set vpage off'
'set grads off'
'set parea 0.0 11.0 0.0 8.5'
'set background 37'
'c'
'set xlab off'
'set ylab off'
'set rbcols 20 21 22 23 24 25 26 27 28 29 31 32 33 34 35 36 37 38 39 40
41 42 43 44 45 46 47 48 49 50 51 52 '
* loop on time steps where each loop does
'set gxout stat'
'd HTSGWsfc'
'set grads off'
'set rbcols 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 20 21
22 23 24 25 26 27 28 29 31 32 33 34 35 36'
'set clevs .5 1 1.5 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 7.5 8 8.5 9
9.5 10 11 12 13 14 15 16 17 18 19 20'
'set gxout barb'
'd (HTSGWsfc*3)'
'run basemap L 3 0 hi '
'set ccolor 0'
'set digsiz .05'
'set digsiz .01'
'd skip(UGRDsfc*1.9, 0);skip(VGRDsfc*1.9, 0)'
'printim 'height_image' png x'xpix' y'ypix
* --- Plot finished!!
* -------------------------------------- 00enp.ctl ----------------
dset /maps/reefcast/wavebin-data/00enp.grb
index /maps/reefcast/wavebin-data/00enp.grb.idx
undef 9.999E+20
title /maps/reefcast/wavebin-data/00enp.grb
* produced by g2ctl v0.0.3b
* griddef=1:0:(373 x 224):grid_template=0:winds(N/S): lat-lon grid:(373
x 224) units 1e-06 input WE:NS output WE:SN res 48 lat 60.500000 to
4.750000 by 0.250000 lon 189.750000 to 282.750000 by 0.250000 #points=83552
dtype grib2
ydef 224 linear 4.750000 0.25
xdef 373 linear 189.750000 0.250000
tdef 61 linear 12Z29dec2020 3hr
zdef 1 linear 1 1
vars 15
DIRPWsfc 0,1,0 10,0,10 ** surface none Primary Wave Direction [deg]
DIRSWsfc 0,1,0 10,0,12 ** surface none Secondary Wave Direction [deg]
HTSGWsfc 0,1,0 10,0,3 ** surface none Significant Height of Combined
Wind Waves and Swell [m]
PERPWsfc 0,1,0 10,0,11 ** surface none Primary Wave Mean Period [s]
PERSWsfc 0,1,0 10,0,13 ** surface none Secondary Wave Mean Period [s]
SWDIR_neg0 0,241,0 10,0,7 ** 0 in sequence none Direction of Swell
Waves [deg]
SWELL_neg0 0,241,0 10,0,8 ** 0 in sequence none Significant Height of
Swell Waves [m]
SWPER_neg0 0,241,0 10,0,9 ** 0 in sequence none Mean Period of Swell
Waves [s]
UGRDsfc 0,1,0 0,2,2 ** surface none U-Component of Wind [m/s]
VGRDsfc 0,1,0 0,2,3 ** surface none V-Component of Wind [m/s]
WDIRsfc 0,1,0 0,2,0 ** surface none Wind Direction (from which
blowing) [deg]
WINDsfc 0,1,0 0,2,1 ** surface none Wind Speed [m/s]
WVDIRsfc 0,1,0 10,0,4 ** surface none Direction of Wind Waves [deg]
WVHGTsfc 0,1,0 10,0,5 ** surface none Significant Height of Wind
Waves [m]
WVPERsfc 0,1,0 10,0,6 ** surface none Mean Period of Wind Waves [s]
ENDVARS
On 12/29/2020 5:25 PM, Stephen Woodbridge wrote:
> Hi All,
>
> I'm generating PNG files in North America and the farther west I go
> the more squashed the image is along the X-axis (longitude) with
> background color to the right and left of the content that gets rendered.
>
>
>
> # East coast:
> ftp://ftpprd.ncep.noaa.gov/pub/data/nccf/com/wave/prod/multi_1.20201229/wna.t00z.grib.grib2
>
>
> # West coast:
> ftp://ftpprd.ncep.noaa.gov/pub/data/nccf/com/wave/prod/multi_1.20201229/enp.t00z.grib.grib2
>
>
> # processing for east coast:
> wgrib2 wna.t00z.grib.grib2 -set_grib_type c3 -grib_out 00wna.grb
>
> g2ctl.pl 00wna.nc > 00wna.ctl
> gribmap -i 00wna.ctl
>
> wgrib2 00.grb -netcdf 00.nc
>
> # loop over 00wna.nc and for each hour generate a PNG
> grads -blc 'run wavecast_new.gs wna '
>
> lat1=20.3
> lat2=48.0
> lon1=262 (-98)
> lon2=302 (-58)
> xpix=8000
> ypix=5540
>
> And I create a world file for each image like:
>
> 0.005
> 0.0
> 0.0
> -0.005
> -98.0
> 48.0
>
>
> Same processing as above but with the enp file for the west coast
> and like wise the western insert has params like:
>
> lat1=30
> lat2=50
> lon1=230 (-130)
> lon2=245 (-115)
> xpix=3000
> ypix=4000
>
> and a world file like:
>
> 0.005
> 0.0
> 0.0
> -0.005
> -130.0
> 50.0
>
> And I get an images like the attached that are overlaid on my map.
>
> I'm thinking that this might be caused by the fact that my image area
> does not align exactly with the area of the grib file. I thought the
> process was extracting data as a sub-windows, but I didn't develop
> these scripts so I'm not sure. If the generated png needs to match the
> grib exactly, how do I extract the dimensions of the grib.
>
> Thanks,
> -Steve
>
>
>
> ga-> q config
> Config: v2.2.1 little-endian readline grib2 netcdf hdf4-sds hdf5
> geotiff shapefile
> Grid Analysis and Display System (GrADS) Version 2.2.1
> Copyright (C) 1988-2018 by George Mason University
> GrADS comes with ABSOLUTELY NO WARRANTY
> See file COPYRIGHT for more information
>
> Configured on 12/29/19 for x86_64-unknown-linux-gnu
>
> This build of GrADS has the following features:
> -+- Byte order is LITTLE ENDIAN
> -+- Athena Widget GUI DISABLED
> -+- Command line editing ENABLED
> -+- GRIB2 interface ENABLED g2clib-1.6.0
> -+- NetCDF interface ENABLED netcdf-4.6.0
> -+- OPeNDAP gridded data interface DISABLED
> -+- OPeNDAP station data interface DISABLED
> -+- HDF4 interface ENABLED hdf-4.2r13
> -+- HDF5 interface ENABLED hdf5-1.10.0
> -+- KML contour output ENABLED
> -+- GeoTIFF and KML/TIFF output ENABLED
> -+- Shapefile interface ENABLED
> The 'q gxconfig' command returns Graphics configuration information
>
>
> $ ncdump -h /maps/reefcast/wavebin-data/00.nc
> netcdf \00 {
> dimensions:
> latitude = 203 ; <<<<<<<< this looks weird for latitude
> longitude = 275 ;
> time = UNLIMITED ; // (61 currently)
> variables:
> double latitude(latitude) ;
> latitude:units = "degrees_north" ;
> latitude:long_name = "latitude" ;
> double longitude(longitude) ;
> longitude:units = "degrees_east" ;
> longitude:long_name = "longitude" ;
> double time(time) ;
> time:units = "seconds since 1970-01-01 00:00:00.0 0:00" ;
> time:long_name = "verification time generated by
> wgrib2 function verftime()" ;
> time:reference_time = 1609243200. ;
> time:reference_time_type = 3 ;
> time:reference_date = "2020.12.29 12:00:00 UTC" ;
> time:reference_time_description = "forecasts or
> accumulated (including analyses), reference date is fixed" ;
> time:time_step_setting = "auto" ;
> time:time_step = 10800. ;
> float WIND_surface(time, latitude, longitude) ;
> WIND_surface:_FillValue = 9.999e+20f ;
> WIND_surface:short_name = "WIND_surface" ;
> WIND_surface:long_name = "Wind Speed" ;
> WIND_surface:level = "surface" ;
> WIND_surface:units = "m/s" ;
> float WDIR_surface(time, latitude, longitude) ;
> WDIR_surface:_FillValue = 9.999e+20f ;
> WDIR_surface:short_name = "WDIR_surface" ;
> WDIR_surface:long_name = "Wind Direction (from which
> blowing)" ;
> WDIR_surface:level = "surface" ;
> WDIR_surface:units = "deg" ;
> float UGRD_surface(time, latitude, longitude) ;
> UGRD_surface:_FillValue = 9.999e+20f ;
> UGRD_surface:short_name = "UGRD_surface" ;
> UGRD_surface:long_name = "U-Component of Wind" ;
> UGRD_surface:level = "surface" ;
> UGRD_surface:units = "m/s" ;
> float VGRD_surface(time, latitude, longitude) ;
> VGRD_surface:_FillValue = 9.999e+20f ;
> VGRD_surface:short_name = "VGRD_surface" ;
> VGRD_surface:long_name = "V-Component of Wind" ;
> VGRD_surface:level = "surface" ;
> VGRD_surface:units = "m/s" ;
> float HTSGW_surface(time, latitude, longitude) ;
> HTSGW_surface:_FillValue = 9.999e+20f ;
> HTSGW_surface:short_name = "HTSGW_surface" ;
> HTSGW_surface:long_name = "Significant Height of
> Combined Wind Waves and Swell" ;
> HTSGW_surface:level = "surface" ;
> HTSGW_surface:units = "m" ;
> float WVPER_surface(time, latitude, longitude) ;
> WVPER_surface:_FillValue = 9.999e+20f ;
> WVPER_surface:short_name = "WVPER_surface" ;
> WVPER_surface:long_name = "Mean Period of Wind Waves" ;
> WVPER_surface:level = "surface" ;
> WVPER_surface:units = "s" ;
> float PERPW_surface(time, latitude, longitude) ;
> PERPW_surface:_FillValue = 9.999e+20f ;
> PERPW_surface:short_name = "PERPW_surface" ;
> PERPW_surface:long_name = "Primary Wave Mean Period" ;
> PERPW_surface:level = "surface" ;
> PERPW_surface:units = "s" ;
> float WVDIR_surface(time, latitude, longitude) ;
> WVDIR_surface:_FillValue = 9.999e+20f ;
> WVDIR_surface:short_name = "WVDIR_surface" ;
> WVDIR_surface:long_name = "Direction of Wind Waves" ;
> WVDIR_surface:level = "surface" ;
> WVDIR_surface:units = "deg" ;
> float DIRPW_surface(time, latitude, longitude) ;
> DIRPW_surface:_FillValue = 9.999e+20f ;
> DIRPW_surface:short_name = "DIRPW_surface" ;
> DIRPW_surface:long_name = "Primary Wave Direction" ;
> DIRPW_surface:level = "surface" ;
> DIRPW_surface:units = "deg" ;
> float WVHGT_surface(time, latitude, longitude) ;
> WVHGT_surface:_FillValue = 9.999e+20f ;
> WVHGT_surface:short_name = "WVHGT_surface" ;
> WVHGT_surface:long_name = "Significant Height of Wind
> Waves" ;
> WVHGT_surface:level = "surface" ;
> WVHGT_surface:units = "m" ;
> float SWELL_0insequence(time, latitude, longitude) ;
> SWELL_0insequence:_FillValue = 9.999e+20f ;
> SWELL_0insequence:short_name = "SWELL_0insequence" ;
> SWELL_0insequence:long_name = "Significant Height of
> Swell Waves" ;
> SWELL_0insequence:level = "0 in sequence" ;
> SWELL_0insequence:units = "m" ;
> float PERSW_surface(time, latitude, longitude) ;
> PERSW_surface:_FillValue = 9.999e+20f ;
> PERSW_surface:short_name = "PERSW_surface" ;
> PERSW_surface:long_name = "Secondary Wave Mean Period" ;
> PERSW_surface:level = "surface" ;
> PERSW_surface:units = "s" ;
> float SWPER_0insequence(time, latitude, longitude) ;
> SWPER_0insequence:_FillValue = 9.999e+20f ;
> SWPER_0insequence:short_name = "SWPER_0insequence" ;
> SWPER_0insequence:long_name = "Mean Period of Swell
> Waves" ;
> SWPER_0insequence:level = "0 in sequence" ;
> SWPER_0insequence:units = "s" ;
> float DIRSW_surface(time, latitude, longitude) ;
> DIRSW_surface:_FillValue = 9.999e+20f ;
> DIRSW_surface:short_name = "DIRSW_surface" ;
> DIRSW_surface:long_name = "Secondary Wave Direction" ;
> DIRSW_surface:level = "surface" ;
> DIRSW_surface:units = "deg" ;
> float SWDIR_0insequence(time, latitude, longitude) ;
> SWDIR_0insequence:_FillValue = 9.999e+20f ;
> SWDIR_0insequence:short_name = "SWDIR_0insequence" ;
> SWDIR_0insequence:long_name = "Direction of Swell
> Waves" ;
> SWDIR_0insequence:level = "0 in sequence" ;
> SWDIR_0insequence:units = "deg" ;
>
> // global attributes:
> :Conventions = "COARDS" ;
> :History = "created by wgrib2" ;
> :GRIB2_grid_template = 0 ;
> }
>
More information about the gradsusr
mailing list