Bug: Extents in PRISM and nldas2 imports are super weird (but daymet extent is spot-on)

[rburghol]

@mwdunlap2004 This of interest as it speaks to the out of memory issue that we get when we don't clip before resampling.

Query to show weird and good extents:

• Good:
  • daymet - select tid, st_astext(st_envelope(rast)) from dh_timeseries_weather where featureid = 617850 and varid = 1454 and extract(year from to_timestamp(tstime)) = 2001 limit 5;
    • POLYGON((-83.66607521306385 44.095589132351506,-73.52855441054511 44.095589132351506,-73.52855441054511 33.95806832983277,-83.66607521306385 33.95806832983277,-83.66607521306385 44.095589132351506))
• Bad:
  • PRISM: select tid, st_astext(st_envelope(rast)) from dh_timeseries_weather where featureid = 617850 and varid = 1453 and extract(year from to_timestamp(tstime)) = 2001 limit 5;
    • POLYGON((-83.64583333330026 44.062499999995055,-41.97916666660025 44.062499999995055,-41.97916666660025 2.395833333295052,-83.64583333330026 2.395833333295052,-83.64583333330026 44.062499999995055))
  • NLDAS2: select tid, st_astext(st_envelope(rast)) from dh_timeseries_weather where featureid = 617850 and varid = 1451 and extract(year from to_timestamp(tstime)) = 2001 limit 5;
    • POLYGON((-83.6255 44.0005,41.3745 44.0005,41.3745 -80.9995,-83.6255 -80.9995,-83.6255 44.0005))
• Code that sets extents (which seems to be correct, and seemed to work for daymet): https://github.com/HARPgroup/meta_model/blob/main/models/raster_met/met/process/02_clip
  • gdalwarp -of "gtiff" -cutline $maskExtent -crop_to_cutline "${originalFile}.repro.gtiff" "${originalFile}${final_ext}"

Short term fix

• I ran these queries to fix these to enable us to avoid query clipping for now, but it will be good when we can try to re-import some PRISM/nldas2 data and see if there was just some problem during original import that has since been fixed.
• prism: update dh_timeseries_weather set rast = st_clip(dh_timeseries_weather.rast, st_envelope(tpr.rast)) from raster_templates as tpr where tpr.varkey = 'daymet_mod_daily' and varid = 1453 ;
  • Note: this query updated 15866 records in around 1-2 minutes total time (I forgot to \timing on unfortunately so I don't have an exact figure).
  • Result of above query after update: POLYGON((-83.64583333330026 44.062499999995055,-73.52083333329216 44.062499999995055,-73.52083333329216 33.937499999986954,-83.64583333330026 33.937499999986954,-83.64583333330026 44.062499999995055))
• nldas2:
  • OriginaL update dh_timeseries_weather set rast = st_clip(dh_timeseries_weather.rast, st_envelope(tpr.rast)) from raster_templates as tpr where tpr.varkey = 'daymet_mod_daily' and varid = 1451 ;
  • modified (update only needed): update dh_timeseries_weather set rast = st_clip(dh_timeseries_weather.rast, st_envelope(tpr.rast)) from raster_templates as tpr where tpr.varkey = 'daymet_mod_daily' and varid = 1451 and NOT (bbox = st_envelope(tpr.rast)) ;
  • UPDATE 385704
  • Time: 681315.766 ms (11:21.316)
  • Result of above query after update: POLYGON((-83.6255 44.0005,-73.5005 44.0005,-73.5005 34.0005,-83.6255 34.0005,-83.6255 44.0005))

Query 1: Update non-conforming bbox:

#!/bin/bash

echo "update dh_timeseries_weather set rast = st_clip(dh_timeseries_weather.rast, st_envelope(tpr.rast)), bbox = st_envelope(tpr.rast) from raster_templates as tpr where tpr.varkey = 'daymet_mod_daily' and varid = 1451 and NOT (bbox = st_envelope(tpr.rast)) ;" | psql -h dbase2 drupal.dh03

[mwdunlap2004]

That's really interesting, I do wonder if just re-importing could fix it, but it's good to know your "short term fix" has led to PRISM being to the correct extent

[rburghol]

Yeah @mwdunlap2004 I'd like to have taken the time to do a reimport but wanted to move forward. Also, the NLDAS2 data is also no fixed, so, we should be able to avoid those clips in future calculations.

[COBrogan]

@rburghol do you want me to look into this a bit? I'm wondering if specifying the source and destination CRS may help. There's a couple of SO posts about vaguely similar issues but in each case they had coordinate transformations not being specified. This _shouldn't__ be a problem in our case, but it's worth a try imo.

It's just weird because that gdalwarp line was called (theoretically) for every data set