search

alexander-petkov / wfas

A placeholder for the WFAS project.
5 stars 1 forks source link

NetCDF data package export #24

Open alexander-petkov opened 4 years ago

alexander-petkov commented 4 years ago

Create a 24-hour data 'package' with temp, RH, solar rad, rainfall,wind speed and direction.

Initial experiments indicate that a 24-hout bundle would be quite large. Would hourly bundles be acceptable?

wmjolly commented 4 years ago

Maybe I should just work directly from the stored data. Alternatively, we could 'tile' the data for smaller files and more efficient processing.

On Thu, Mar 19, 2020 at 5:24 PM alexander-petkov notifications@github.com wrote:

Create a 24-hour data 'package' with temp, RH, solar rad, rainfall,wind speed and direction.

Initial experiments indicate that a 24-hout bundle would be quite large. Would hourly bundles be acceptable?

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/alexander-petkov/wfas/issues/24, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4G3D7JQJSQQP6TI6TTUSLRIKSTLANCNFSM4LPY2EYQ .

alexander-petkov commented 4 years ago

1 hour Netcdf archive is 80MB.

alexander-petkov commented 4 years ago

NetCDF file from windninja output

tail -n +7 windninja_solar.asc|tac |cdo -f nc4 -b I16 \ 
   -settaxis,2019-09-10,22:00:00,hours -setname,solar \
   -input,cephfs/tmp/varanl_grid outfile.nc

The tail command prints the file from the 7 line on, tac reverses it (bacause for some reason cdo flips the result upside down), then cdo uses varanl_grid to set projection in orig Grib files, sets the variable name, the time stamp and period, output data type, and file format to NetCDF4

alexander-petkov commented 4 years ago

Here is information from 1-hour RTMA archive with temp, RH, solar rad, rainfall,wind speed and direction

cdo sinfon out.nc 
   File format : NetCDF4
    -1 : Institut Source   T Steptype Levels Num    Points Num Dtype : Parameter name
     1 : unknown  unknown  v instant       1   1   3744965   1  F32  : 2t            
     2 : unknown  unknown  v instant       1   2   3744965   1  F32  : 10wdir        
     3 : unknown  unknown  v instant       1   2   3744965   1  F32  : 10si          
     4 : unknown  unknown  v instant       1   1   3744965   1  F32  : 2r            
     5 : unknown  unknown  v instant       1   3   2953665   2  F32  : tp            
     6 : unknown  unknown  v instant       1   3   3744965   1  I16  : solar         
   Grid coordinates :
     1 : projection               : points=3744965 (2345x1597)
                          mapping : lambert_conformal_conic
                                x : 0 to 5953064 by 2539.703 m
                                y : 0 to 4053366 by 2539.703 m
     2 : projection               : points=2953665 (2145x1377)
                          mapping : lambert_conformal_conic
                              x_2 : 0 to 5445123 by 2539.703 m
                              y_2 : 0 to 3494631 by 2539.703 m
   Vertical coordinates :
     1 : height                   : levels=1
                           height : 2 m
     2 : height                   : levels=1
                         height_2 : 10 m
     3 : surface                  : levels=1
   Time coordinate :  1 step
     RefTime =  2019-09-10 22:00:00  Units = hours  Calendar = proleptic_gregorian
  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss
  2019-09-10 22:00:00
cdo sinfon: Processed 6 variables over 1 timestep [0.01s 53MB]

Total precip (tp) is not on the same grid, as the rest of the variables.
Would that be acceptable?

EDIT: It is possible to regrid Precip data to the same grid as the rest of the variables. However, it takes significantly longer to assemble hourly archive:

time cdo -P 4 -f nc4  -merge -selname,2t,10wdir,10si cephfs/wfas/data/rtma/varanl/grb/rtma2p5.20190910/rtma2p5.t22z.2dvaranl_ndfd.grb2_wexp -setgrid,cephfs/tmp/varanl_grid cephfs/wfas/data/rtma/rhm/rtma2p5.20190910/rtma2p5.t22z.2dvaranl_ndfd.grb2_wexp -remapycon,cephfs/tmp/varanl_grid cephfs/wfas/data/rtma/pcp/rtma2p5.20190910/rtma2p5.2019091022.pcp.184.grb2 outfile.nc out.nc
cdo(2) selname: Process started
cdo(3) setgrid: Process started
cdo(4) remapcon: Process started
cdo(4) remapcon: YAC first order conservative weights from curvilinear (2145x1377) to projection (2345x1597) grid, with source mask (1558737)
cdo(2) selname: Processed 11234895 values from 13 variables over 1 timestep
cdo(3) setgrid: Processed 3744965 values from 1 variable over 1 timestep
cdo(4) remapcon: Processed 2953665 values from 1 variable over 1 timestep
cdo merge: Processed 22469790 values from 6 variables over 1 timestep [82.81s 2053MB]

real    1m22.927s
user 4m20.375s
sys 0m1.823s

Example with all variables on the same grid:

cdo sinfon out.nc 
   File format : NetCDF4
    -1 : Institut Source   T Steptype Levels Num    Points Num Dtype : Parameter name
     1 : unknown  unknown  v instant       1   1   3744965   1  F32  : 2t            
     2 : unknown  unknown  v instant       1   2   3744965   1  F32  : 10wdir        
     3 : unknown  unknown  v instant       1   2   3744965   1  F32  : 10si          
     4 : unknown  unknown  v instant       1   1   3744965   1  F32  : 2r            
     5 : unknown  unknown  v instant       1   3   3744965   1  F32  : tp            
     6 : unknown  unknown  v instant       1   3   3744965   1  I16  : solar         
   Grid coordinates :
     1 : projection               : points=3744965 (2345x1597)
                          mapping : lambert_conformal_conic
                                x : 0 to 5953064 by 2539.703 m
                                y : 0 to 4053366 by 2539.703 m
   Vertical coordinates :
     1 : height                   : levels=1
                           height : 2 m
     2 : height                   : levels=1
                         height_2 : 10 m
     3 : surface                  : levels=1
   Time coordinate :  1 step
     RefTime =  2019-09-10 22:00:00  Units = hours  Calendar = proleptic_gregorian
  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss
  2019-09-10 22:00:00
cdo sinfon: Processed 6 variables over 1 timestep [0.02s 53MB]

Without regridding it takes 8 seconds.

wmjolly commented 4 years ago

Can we resample the precip to match?

On Fri, Mar 20, 2020, 16:06 alexander-petkov notifications@github.com wrote:

Here is information from 1-hour RTMA archive with temp, RH, solar rad, rainfall,wind speed and direction

cdo sinfon out.nc File format : NetCDF4 -1 : Institut Source T Steptype Levels Num Points Num Dtype : Parameter name 1 : unknown unknown v instant 1 1 3744965 1 F32 : 2t 2 : unknown unknown v instant 1 2 3744965 1 F32 : 10wdir 3 : unknown unknown v instant 1 2 3744965 1 F32 : 10si 4 : unknown unknown v instant 1 1 3744965 1 F32 : 2r 5 : unknown unknown v instant 1 3 2953665 2 F32 : tp 6 : unknown unknown v instant 1 3 3744965 1 I16 : solar Grid coordinates : 1 : projection : points=3744965 (2345x1597) mapping : lambert_conformal_conic x : 0 to 5953064 by 2539.703 m y : 0 to 4053366 by 2539.703 m 2 : projection : points=2953665 (2145x1377) mapping : lambert_conformal_conic x_2 : 0 to 5445123 by 2539.703 m y_2 : 0 to 3494631 by 2539.703 m Vertical coordinates : 1 : height : levels=1 height : 2 m 2 : height : levels=1 height_2 : 10 m 3 : surface : levels=1 Time coordinate : 1 step RefTime = 2019-09-10 22:00:00 Units = hours Calendar = proleptic_gregorian YYYY-MM-DD hh:mm:ss YYYY-MM-DD hh:mm:ss YYYY-MM-DD hh:mm:ss YYYY-MM-DD hh:mm:ss 2019-09-10 22:00:00 cdo sinfon: Processed 6 variables over 1 timestep [0.01s 53MB]

Total precip (tp) is not on the same grid, as the rest of the variables. Would that be acceptable?

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/alexander-petkov/wfas/issues/24#issuecomment-601930248, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4G3D3IDIPOQJBDA6TZBRLRIPSEXANCNFSM4LPY2EYQ .

alexander-petkov commented 4 years ago

Can we resample the precip to match?

Yes--putting Precip on the same grid takes longer, about minute and a half for each time step. But it's nice and uniform.

alexander-petkov commented 4 years ago

This NetCDF archive is for 16 hours, about 1.4G file size Screenshot from 2020-03-20 23-17-39

A full 24 hour archive should be about 2G. The slowest part of the process is the regridding of Total Precip.

alexander-petkov commented 4 years ago

This is now working for RTMA.

NDFD doesn't have continuous 24-hour coverage. There are gaps in the data of 2, 3 or more hours.

wmjolly commented 4 years ago

I know the files will be big but 24 hours of data will make processing overall much faster and easier.

On Fri, Mar 20, 2020 at 11:22 PM alexander-petkov notifications@github.com wrote:

This NetCDF archive is for 16 hours, about 1.4G file size [image: Screenshot from 2020-03-20 23-17-39] https://user-images.githubusercontent.com/39599557/77220027-539a9380-6b01-11ea-8065-5679850a16b6.png

A full 24 hour archive should be about 2G. The slowest part of the process is the regridding of Total Precip.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/alexander-petkov/wfas/issues/24#issuecomment-601997131, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4G3D2WTS5Z7A2GO7FDCK3RIRFINANCNFSM4LPY2EYQ .

wmjolly commented 4 years ago

So, the NetCDF files have the "height" field in there 24 times. I think we can extract that field once, store in the database and not include it in the NetCDF files in the future.

On Mon, Mar 23, 2020 at 9:54 AM Matt Jolly wmjolly@gmail.com wrote:

I know the files will be big but 24 hours of data will make processing overall much faster and easier.

On Fri, Mar 20, 2020 at 11:22 PM alexander-petkov < notifications@github.com> wrote:

This NetCDF archive is for 16 hours, about 1.4G file size [image: Screenshot from 2020-03-20 23-17-39] https://user-images.githubusercontent.com/39599557/77220027-539a9380-6b01-11ea-8065-5679850a16b6.png

A full 24 hour archive should be about 2G. The slowest part of the process is the regridding of Total Precip.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/alexander-petkov/wfas/issues/24#issuecomment-601997131, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4G3D2WTS5Z7A2GO7FDCK3RIRFINANCNFSM4LPY2EYQ .

alexander-petkov commented 4 years ago

So, the NetCDF files have the "height" field in there 24 times. I think we can extract that field once, store in the database and not include it in the NetCDF files in the future.

What are the steps to reproduce this behaviour?

wmjolly commented 4 years ago

There is a height and height_2 value. Are they just 2-d arrays and not 3-d like the other vars?

wmjolly@wmjolly-VirtualBox:~/Downloads$ ncdump -h rtma2p5.20200320.nc netcdf rtma2p5.20200320 { dimensions: time = UNLIMITED ; // (24 currently) x = 2345 ; y = 1597 ; height = 1 ; height_2 = 1 ; variables: double time(time) ; time:standard_name = "time" ; time:units = "hours since 2020-3-20 00:00:00" ; time:calendar = "proleptic_gregorian" ; time:axis = "T" ; double x(x) ; x:standard_name = "projection_x_coordinate" ; x:units = "m" ; x:axis = "X" ; double y(y) ; y:standard_name = "projection_y_coordinate" ; y:units = "m" ; y:axis = "Y" ; int Lambert_Conformal ; Lambert_Conformal:grid_mapping_name = "lambert_conformal_conic" ; Lambert_Conformal:standard_parallel = 25. ; Lambert_Conformal:longitude_of_central_meridian = 265. ; Lambert_Conformal:latitude_of_projection_origin = 25. ; Lambert_Conformal:earth_radius = 6367470. ; Lambert_Conformal:false_easting = 3269232.20488998 ; Lambert_Conformal:false_northing = 263643.063253065 ; Lambert_Conformal:longitudeOfFirstGridPointInDegrees = 233.723448 ; Lambert_Conformal:latitudeOfFirstGridPointInDegrees = 19.228976 ; double height(height) ; height:standard_name = "height" ; height:long_name = "height" ; height:units = "m" ; height:positive = "up" ; height:axis = "Z" ; double height_2(height_2) ; height_2:standard_name = "height" ; height_2:long_name = "height" ; height_2:units = "m" ; height_2:positive = "up" ; height_2:axis = "Z" ; float \2t(time, height, y, x) ; \2t:standard_name = "air_temperature" ; \2t:long_name = "2 metre temperature" ; \2t:units = "K" ; \2t:param = "0.0.0" ; \2t:grid_mapping = "Lambert_Conformal" ; float \10wdir(time, height_2, y, x) ; \10wdir:long_name = "10 metre wind direction" ; \10wdir:units = "Degree true" ; \10wdir:param = "0.2.0" ; \10wdir:grid_mapping = "Lambert_Conformal" ; float \10si(time, height_2, y, x) ; \10si:long_name = "10 metre wind speed" ; \10si:units = "m s-1" ; \10si:param = "1.2.0" ; \10si:grid_mapping = "Lambert_Conformal" ; float \2r(time, height, y, x) ; \2r:standard_name = "relative_humidity" ; \2r:long_name = "2 metre relative humidity" ; \2r:units = "%" ; \2r:param = "1.1.0" ; \2r:grid_mapping = "Lambert_Conformal" ; float tp(time, y, x) ; tp:long_name = "Total Precipitation" ; tp:units = "kg m-2" ; tp:param = "8.1.0" ; tp:grid_mapping = "Lambert_Conformal" ; tp:_FillValue = -9.e+33f ; tp:missing_value = -9.e+33f ; float solar(time, y, x) ; solar:grid_mapping = "Lambert_Conformal" ;

// global attributes: :CDI = "Climate Data Interface version 1.9.7.1 (http://mpimet.mpg.de/cdi)" ; :history = "Mon Mar 23 01:12:18 2020: cdo -P 4 mergetime /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t00z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t01z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t02z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t03z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t04z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t05z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t06z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t07z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t08z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t09z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t10z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t11z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t12z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t13z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t14z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t15z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t16z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t17z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t18z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t19z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t20z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t21z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t22z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t23z.2dvaranl_ndfd.grb2_wexp.nc /mnt/cephfs/netcdf/rtma/ rtma2p5.20200320.nc\nMon Mar 23 00:29:27 2020: cdo -P 4 -f nc4 -merge -selname,2t,10wdir,10si /mnt/cephfs/wfas/data/rtma/varanl/grb/rtma2p5.20200320/rtma2p5.t00z.2dvaranl_ndfd.grb2_wexp -setgrid,/mnt/cephfs/wfas/data/rtma/varanl/varanl_grid /mnt/cephfs/wfas/data/rtma/rhm/grb/rtma2p5.20200320/rtma2p5.t00z.2dvaranl_ndfd.grb2_wexp -remapycon,/mnt/cephfs/wfas/data/rtma/varanl/varanl_grid /mnt/cephfs/wfas/data/rtma/pcp/grb/rtma2p5.20200320/rtma2p5.2020032000.pcp.184.grb2 /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/rtma2p5.t00z.2dvaranl_ndfd.grb2_wexp.solar /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/ rtma2p5.t00z.2dvaranl_ndfd.grb2_wexp.nc\nMon Mar 23 00:29:25 2020: cdo -P 4 -f nc4 -settaxis,2020-03-20,00:00:00,hours -setname,solar -input,/mnt/cephfs/wfas/data/rtma/varanl/varanl_grid /mnt/cephfs/netcdf/rtma/rtma2p5.20200320/rtma2p5.t00z.2dvaranl_ndfd.grb2_wexp.solar" ; :_NCProperties = "version=2,netcdf=4.6.2,hdf5=1.10.5" ; :institution = "ECMWF" ; :Conventions = "CF-1.6" ; :cdo_openmp_thread_number = 4 ; :CDO = "Climate Data Operators version 1.9.7.1 (http://mpimet.mpg.de/cdo)" ; }

On Mon, Mar 30, 2020 at 10:13 AM alexander-petkov notifications@github.com wrote:

So, the NetCDF files have the "height" field in there 24 times. I think we can extract that field once, store in the database and not include it in the NetCDF files in the future.

What are the steps to reproduce this behaviour?

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/alexander-petkov/wfas/issues/24#issuecomment-606095604, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4G3D6726AKBZS2MQJN4CLRKDALLANCNFSM4LPY2EYQ .

alexander-petkov commented 4 years ago

There is a height and height_2 value. Are they just 2-d arrays and not 3-d like the other vars?

Each height dimension should be a single value--2 and 10 m respectively. This from ncdump:

height = 2 ;

height_2 = 10 ;

alexander-petkov commented 4 years ago

There is an issue with using GDAL to read NetCDF files with 2 or more bands.

If I create a single band NetCDF file, CRS information is displayed correctly. Adding another band results in error:

cdo -P 4 -f nc4  -selname,2t /mnt/cephfs/wfas/data/rtma/varanl/grb/rtma2p5.20200324/rtma2p5.t00z.2dvaranl_ndfd.grb2_wexp 2t.nc

gdalsrsinfo 2t.nc 

PROJ.4 : +proj=lcc +lat_1=25 +lat_0=25 +lon_0=265 +k_0=1 +x_0=3269232.20488998 +y_0=263643.063253065 +R=6367470 +units=m +no_defs

cdo -P 4 -f nc4  -selname,2t,10wdir /mnt/cephfs/wfas/data/rtma/varanl/grb/rtma2p5.20200324/rtma2p5.t00z.2dvaranl_ndfd.grb2_wexp 2t.nc

gdalsrsinfo 2t.nc 
ERROR 1: ERROR - failed to load SRS definition from 2t.nc

This is apparently a documented issue: https://issues.qgis.org/issues/21783

wmjolly commented 4 years ago

Copy. My mistake. I interpreted the ncdump -h as showing a layer each timestep.

MJ

On Mon, Mar 30, 2020 at 10:48 AM alexander-petkov notifications@github.com wrote:

There is a height and height_2 value. Are they just 2-d arrays and not 3-d like the other vars?

Each height dimension should be a single value--2 and 10 m respectively. This from ncdump:

height = 2 ;

height_2 = 10 ;

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/alexander-petkov/wfas/issues/24#issuecomment-606114250, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4G3D5QKDOT5AGZMUVVXWLRKDEOJANCNFSM4LPY2EYQ .

alexander-petkov commented 4 years ago

GDAL apparently reads Spatial reference info on a per-dataset basis :

gdalsrsinfo NETCDF:"/mnt/cephfs/netcdf/rtma/out.nc":10wdir

PROJ.4 : +proj=lcc +lat_1=25 +lat_0=25 +lon_0=265 +k_0=1 +x_0=3269232.20488998 +y_0=263643.063253065 +R=6367470 +units=m +no_defs

OGC WKT2:2018 :
PROJCRS["unnamed",
    BASEGEOGCRS["unknown",
        DATUM["unnamed",
            ELLIPSOID["Sphere",6367470,0,
                LENGTHUNIT["metre",1,
                    ID["EPSG",9001]]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]],
    CONVERSION["unnamed",
        METHOD["Lambert Conic Conformal (1SP)",
            ID["EPSG",9801]],
        PARAMETER["Latitude of natural origin",25,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",265,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",1,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",3269232.20488998,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",263643.063253065,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]],
        PARAMETER["standard_parallel_1",25,
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]],
    CS[Cartesian,2],
        AXIS["easting",east,
            ORDER[1],
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]],
        AXIS["northing",north,
            ORDER[2],
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]]]

Maybe that's how I should load a NetCDF dataset in QGIS?

Edit: Huh, that worked.... Although the data is flipped.

Screen Shot 2020-03-31 at 1 46 01 PM

alexander-petkov commented 4 years ago

Change the cdo command such that over the top history is not included in the metadata.

alexander-petkov commented 4 years ago

I believe I have this buttoned up--all data sets are oriented properly. I have left the original CRS definition as in the RTMA archive, and will change it if there are problems. I didn't notice any in QGIS or Geoserver.

I will rerun the script again to regenerate the NetCDF bundles

Initial commit: https://github.com/alexander-petkov/wfas/commit/ba5cce0d8d5b4991c9bbb93776adbf64248d5ee1

wmjolly commented 4 years ago

Send me the link when you have a new file for me to work with. Thanks Alex.

On Tue, Mar 31, 2020, 15:44 alexander-petkov notifications@github.com wrote:

I believe I have this buttoned up--all data sets are oriented properly. I have left the original CRS definition as in the RTMA archive, and will change it if there are problems. I didn't notice any in QGIS or Geoserver.

I will rerun the script again to regenerate the NetCDF bundles

Initial commit: ba5cce0 https://github.com/alexander-petkov/wfas/commit/ba5cce0d8d5b4991c9bbb93776adbf64248d5ee1

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/alexander-petkov/wfas/issues/24#issuecomment-606893408, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4G3DZ4IYKHODRWTC2UVJDRKJP5JANCNFSM4LPY2EYQ .

alexander-petkov commented 4 years ago

I modifed the grid definition, to align the data in -180 to 180 longitude space:

#
# gridID 1
#
gridtype  = projection
gridsize  = 3744965
xsize     = 2345
ysize     = 1597
xname     = x
xunits    = "m"
yname     = y
yunits    = "m"
xfirst    = -3272417.140
xinc      = 2539.703
yfirst    = -265067.354
yinc      = 2539.703
grid_mapping = Lambert_Conformal
grid_mapping_name = lambert_conformal_conic
standard_parallel = 25.
longitude_of_central_meridian = -95.
latitude_of_projection_origin = 25.
earth_radius = 6367470.
false_easting = 0
false_northing = 0
alexander-petkov commented 3 years ago

I had an epiphany (aka DUH! moment) regarding the NetCDF package export.

It all stemmed from this conversation:

I need a way to easily get the GFS forecast as a NetCDF package.
Same as the RTMA>
Same as the RTMA one you made.
I can request the individual variables over WCS
and get the NetCDF back.
But I'm thinking a 'modeling package' would be easier.
...
Actually, I want them ALL in WGS84 in the future.
It's a pain for modeling because I still need to know the latitude for each grid cell.

So, why not use Geoserver itself for getting NetCDF data??

Pros:

  1. It is possible to use subsetting on any dimension (Geoserver related notes: WCS Requests)
  2. It is possible to request a coverage in a different projection.
  3. It is possible to resample the pixel resolution
  4. Variables missing in some of the original source data (ex Solar Radiation) are already calculated
  5. Accessible from anywhere.

Cons

  1. Original data might be distorted, after converting to a Geotiff format, or reprojecting upon export.
wmjolly commented 3 years ago

Wonder if we could just make a WPS that wraps the multiple WCS requests that packages all the needed variables into a single call?

On Thu, Mar 11, 2021, 01:12 alexander-petkov @.***> wrote:

I had an epiphany (aka DUH! moment) regarding the NetCDF package export.

It all stemmed from this conversation:

I need a way to easily get the GFS forecast as a NetCDF package. Same as the RTMA> Same as the RTMA one you made. I can request the individual variables over WCS and get the NetCDF back. But I'm thinking a 'modeling package' would be easier. ... Actually, I want them ALL in WGS84 in the future. It's a pain for modeling because I still need to know the latitude for each grid cell.

So, why not use Geoserver itself for getting NetCDF data??

Pros:

  1. It is possible to use subsetting on any dimension (Geoserver related notes: WCS Requests https://github.com/alexander-petkov/wfas/wiki/Geoserver-related-notes#wcs-requests )
  2. It is possible to request a coverage in a different projection.
  3. It is possible to resample the pixel resolution
  4. Variables missing in some of the original source data (ex Solar Radiation) is already calculated
  5. Accessible from anywhere.

Cons

  1. Original data might be distorted, after converting to a Geotiff format, or reprojecting upon export.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/alexander-petkov/wfas/issues/24#issuecomment-796552709, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4G3DYXTXGWFMWPMHMU2O3TDBUMTANCNFSM4LPY2EYQ .

alexander-petkov commented 3 years ago

Here is a test of the idea above, using the GFS archive and specifying a 24-hour time period and a spatial subset that covers California:

vars=(Temperature Wind_direction Wind_speed Relative_humidity Total_precipitation Solar_radiation)
for var in ${vars[@]};do wget -q --no-check-certificate "https://aws.wfas.net/geoserver/gfs/wcs?service=WCS&version=2.0.1&request=GetCoverage&coverageId=gfs:${var}&outputCRS=http://www.opengis.net/def/crs/EPSG/0/4326&format=application/x-netcdf&subset=time(%222021-03-22T00:00:00Z%22,%222021-03-22T23:59:59Z%22)&subset=Lat(30,42)&subset=Long(-125,-115)" -O gfs_${var}.nc ;done
ls gfs_*.nc
gfs_Relative_humidity.nc  gfs_Temperature.nc          gfs_Wind_direction.nc
gfs_Solar_radiation.nc    gfs_Total_precipitation.nc  gfs_Wind_speed.nc

cdo -P 4 -O -f nc4 -setgatt,history,"merged all timesteps" -merge gfs_*.nc gfs.nc
cdo(1) merge: Process started
cdo(1) merge: Processed 92160 values from 6 variables over 48 timesteps.
cdo    setgatt: Processed 92160 values from 6 variables over 8 timesteps [0.01s 58MB].

Result: Screenshot from 2021-03-11 01-20-54

Edit: First problem I see is that height dimension is lost on some variables, but this is because coverages were configured without Elevation axis.

alexander-petkov commented 3 years ago

Remap RTMA to a common grid from Geoserver-produced NetCDF files

The first command extracts the grid from rtma_Temperature.nc file, and the second remaps Total precip to a common grid, and merges all files into one bundle:

cdo griddes rtma_Temperature.nc  >rtma.grid
cdo -P 4 -O -f nc4 -setgrid,rtma.grid \
   -setgatt,history,"merged all timesteps" 
   -merge rtma_Temperature.nc rtma_Relative_humidity.nc \
   rtma_Solar_radiation.nc -remapcon,rtma.grid rtma_Total_precipitation.nc \
   rtma_Wind_direction.nc rtma_Wind_speed.nc  rtma.nc

Result: Screenshot from 2021-03-12 01-58-01

alexander-petkov commented 3 years ago

Coordinates and axis names need to be provided in the Coverage's original projection, even if requested in another projection. Example: geoserver/wcs?service=WCS&version=2.0.1&request=GetCoverage&outputCRS=http://www.opengis.net/def/crs/EPSG/0/4326&format=application/x-netcdf&subset=time(%222021-07-20T00:00:00Z%22,%222021-07-20T23:59:59Z%22)&subset=y(0,90838.5)&subset=x(0,83186.5)&coverageId=rtma:Total_precipitation