Freshwater-Initiative / SkagitLandslideHazards

Seattle City Light is interested in improving understanding of landslide hazard and sediment transport to ensure reliable and cost-effective hydropower generation.
4 stars 4 forks source link

DHSVM groundwater output #1

Open ChristinaB opened 6 years ago

ChristinaB commented 6 years ago

The output for the Skagit model is uploaded in the Code

ChristinaB commented 6 years ago
ChristinaB commented 5 years ago

Not using 1915-2010 data - using SeYeun's model time periods instead.

ChristinaB commented 5 years ago
ChristinaB commented 5 years ago

Number of Rows = 1020 # Number of rows Number of Columns = 916 # Number of columns depth to groundwater CNRM has 62 years/maps for depth to groundwater = 63240 rows in the binary file ./myconvert float ascii Map.Soil.TableDepth.bin Map.Soil.TableDepth.asc 63240 916 CNRM45: Error reading input: Row: 59160 Column: 0 - > 58 years only printed? CNRM85: Error reading input: Row: 60180 Column: 0 - > 59 years only printed?

snow CNRM has 61 years/maps for snow = 62220 rows in the binary file ./myconvert float ascii Map.Snow.Swq.bin Map.Snow.Swq.asc 62220 916 CNRM45: Error reading input: Row: 58140 Column: 0 - > 57 years only printed? CNRM85: Error reading input: Row: 58140 Column: 0 - > 57 years only printed?

ice CNRM has 61 years/maps for ice = 62220 rows in the binary file ./myconvert float ascii Map.Snow.Iwq.bin Map.Snow.Iwq.asc 62220 916 CNRM45: Error reading input: Row: 59160 Column: 0 - > 58 years only printed? CNRM85: Error reading input: Row: 58140 Column: 0 - > 57 years only printed?

ChristinaB commented 5 years ago

SBATCH --time=300:00:00

in the SLURM submissin needs to be extended to 400 hours or rerun the last 10 years

ChristinaB commented 5 years ago

Draft matlab code for processing and checking that depth to water table values work as expected. See C:\Users\cband\Skagit\SCLlandsliding\DHSVMoutputs_20190405\output_SCLlandslide_bclivlow_G_CNRM-CM5__rcp45 nrow=1020; [nrowtotal, ncol]=size(watertabledata) mapcount=nrowtotal/nrow;

%Print Livneh data cd C:\Users\cband\Skagit\SCLlandsliding\DHSVMoutputs_20190405\output_SCLlandslide_bclivlow_G_CNRM-CM5__rcp45

watertabledata=load('Map.Soil.TableDepth.asc'); snowdata=load('Map.Snow.Swq.asc'); icedata=load('Map.Snow.Swq.asc');

dateswatertableCNRM45=load('dates-watertableCNRM45.txt');

% WATER TABLE DEPTH check length of data read in i=1 datafilename=[{strcat('watertabledepth',num2str(i),'-',num2str(dateswatertableCNRM45(i)))}]; for i=2:mapcount datafilename=[datafilename {strcat('watertabledepth',num2str(i),'-',num2str(dateswatertableCNRM45(i)))}]; end

GCM='CNRM-CM5__rcp45';

%Print depth to groundwater data %for i=1:2 dlmwrite(strcat(char(GCM),'',char(datafilename(1)),'.asc'),watertabledata(1:nrow,:),'delimiter',' ','precision','%.1f'); dlmwrite(strcat(char(GCM),'',char(datafilename(2)),'.asc'),watertabledata(nrow+1:nrow*2,:),'delimiter',' ','precision','%.1f');

for i=3:mapcount dlmwrite(strcat(char(GCM),'_',char(datafilename(i)),'.asc'),watertabledata(nrow(i-1)+1:nrowi,:),'delimiter',' ','precision','%.1f'); end

% SNOW & ICE check length of data read in i=1 datafilename=[{strcat('snowdata',num2str(i),'-',num2str(dateswatertableCNRM45(i)))}]; datafilename=[{strcat('icedata',num2str(i),'-',num2str(dateswatertableCNRM45(i)))}]; for i=2:mapcount datafilename=[datafilename {strcat('snowdata',num2str(i),'-',num2str(dateswatertableCNRM45(i)))}]; datafilename=[datafilename {strcat('icedata',num2str(i),'-',num2str(dateswatertableCNRM45(i)))}]; end

GCM='CNRM-CM5__rcp45';

%Print depth to groundwater data %for i=1:2 dlmwrite(strcat(char(GCM),'',char(datafilename(1)),'.asc'),watertabledata(1:nrow,:),'delimiter',' ','precision','%.1f'); dlmwrite(strcat(char(GCM),'',char(datafilename(2)),'.asc'),watertabledata(nrow+1:nrow2,:),'delimiter',' ','precision','%.1f'); dlmwrite(strcat(char(GCM),'',char(datafilename(1)),'.asc'),snowdata(1:nrow,:),'delimiter',' ','precision','%.1f'); dlmwrite(strcat(char(GCM),'',char(datafilename(2)),'.asc'),snowdata(nrow+1:nrow2,:),'delimiter',' ','precision','%.1f'); dlmwrite(strcat(char(GCM),'',char(datafilename(1)),'.asc'),icedata(1:nrow,:),'delimiter',' ','precision','%.1f'); dlmwrite(strcat(char(GCM),'',char(datafilename(2)),'.asc'),icedata(nrow+1:nrow*2,:),'delimiter',' ','precision','%.1f');

for i=3:mapcount dlmwrite(strcat(char(GCM),'',char(datafilename(i)),'.asc'),snowdata(nrow(i-1)+1:nrowi,:),'delimiter',' ','precision','%.1f'); dlmwrite(strcat(char(GCM),'',char(datafilename(i)),'.asc'),icedata(nrow(i-1)+1:nrowi,:),'delimiter',' ','precision','%.1f'); end

SaveAsciiRaster(varname, header)

ChristinaB commented 5 years ago

testingoutputs_20190405.xlsx

It works!!

RondaStrauch commented 5 years ago

@ChristinaB What does this mean in English? What worked?

ChristinaB commented 5 years ago

image

It works! the saturation events from this GCM shift from peaks in June to peaks in winter (relatively speaking over 30 year intervals).

ChristinaB commented 5 years ago

To Do by May 2, 2019

Steps: Read in tarball

Read in header

Loop through files

Clip to model extent

Calculate a probability distribution of each time period at each grid cell (150 m).

Attach to landlab raster model grid

Use the same distribution at 30 m resolution

Coupling complete! Test input to landlab landsliding component.

ChristinaB commented 5 years ago

Here is code for processing binary grid outputs into ascii files

Number of Rows = 1020 # Number of rows Number of Columns = 916 # Number of columns

Historic 44 year for depth to groundwater, 43 years for ice and snow ./myconvert float ascii Map.Soil.TableDepth.bin Map.Soil.TableDepth.asc 44880 916 ./myconvert float ascii Map.Snow.Swq.bin Map.Snow.Swq.asc 43860 916 ./myconvert float ascii Map.Snow.Iwq.bin Map.Snow.Iwq.asc 43860 916

depth to groundwater CNRM has 62 years/maps for depth to groundwater = 63240 rows in the binary file ./myconvert float ascii Map.Soil.TableDepth.bin Map.Soil.TableDepth.asc 63240 916 ./myconvert float ascii Map.Snow.Swq.bin Map.Snow.Swq.asc 62220 916 ./myconvert float ascii Map.Snow.Iwq.bin Map.Snow.Iwq.asc 62220 916

ChristinaB commented 5 years ago

@NCristea The saturation extent dates are on Github, and depth to water table DHSVM outputs on HydroShare. Download dtw_DHSVMoutput_20190502.tar

ChristinaB commented 5 years ago

@NCristea Hyak shuts down on May 14, so it won't let me start a two week run that won't finish in time. I am running for the next 150 hours, which may get a good way into the new runs. But we need to just develop the code to work from Map.Depth.X with inputs of rows, columns, and number of maps, and convert into ascii with names by date (see above). The final results won't be ready until end of May.

ChristinaB commented 5 years ago

https://github.com/Freshwater-Initiative/OGH_group_tutorial/blob/master/Observatory_usecase7_xmapLandlab.ipynb

See oxl.calculateUTMbounds() and oxl.mappingfileToRaster() functions in https://github.com/Freshwater-Initiative/Observatory/blob/master/ogh/ogh_xarray_landlab.py

ChristinaB commented 5 years ago

@NCristea you can download a zipped file of raster and ascii files and headers at this link: SCL Domain

Extent Resolution/Grid Cell Size xllcorner yllcorner ncols nrows
Observatory Climate Forcings SCL extent 6 km 591315 5354895 13 13
Skagit Hydrology & Climate Change Study Domain 150 meter 535665 5311395 916 1020
Seattle City Light Skagit Study Domain - Hydrology 150 meter 591315 5354895 545 544
Seattle City Light Skagit Study Domain - Landslide/Fire 30 meter 591315 5354895 2725 2720
ChristinaB commented 4 years ago

https://github.com/Freshwater-Initiative/SkagitLandslideHazards/blob/master/dhsvm_config_files/INPUT.skagitSCLlandslide.historic1960.bclivlow.G

image

Let's use this list to verify the dates in NetCDF DHSVM output/Landlab input using the Skagit input config files and Skagit data on HydroShare in GIS folder. https://www.hydroshare.org/resource/70d746c7da584ae6bd2f88deb5a4c188/

ChristinaB commented 4 years ago

Time 1: dates from Sauk extent model run dtw_DHSVMoutputs_20190405.tar Model Runs Labeled 20190405 match the config files on my computer are dated March 8. This is what is on Github in comment above. image

Time 2: dates from Skagit extent model run 20190502. also correspond. This file is on Hyak. image

Time 3: re-compressing and sending outputs from Hyak processed in May --dtw_DHSVMoutputs_20190925.tar

ChristinaB commented 4 years ago

To Do - Run model for dates from Sauk extent (done - spatially representative for rain and rain-snow dominated area), Skagit extent (will bias toward snow melt days for max saturated event), GoodellCreek extent (will be specific to the fire; rain-snow focus). Is the landslide model sensitive? Does all of this matter? Or is Queen Slope the Master of the Slippery Future Universe?

ChristinaB commented 4 years ago

@NCristea We flipped the Y axis of the xarray netcdf. data = xr.open_dataset('dtw_historic_with_dates.nc') flippy = data.reindex(y=list(reversed(data.y))) this won't plot flipped in xarray plotting because it is too smart and fixes the flipped UTM coordinates in the map. But investigating the sequence shows the correct Y sequence is flipped. The first Y coordinate should be 5354910. We have this working and validated it using the Landlab grid building/netcdf read writing.

image

ChristinaB commented 4 years ago

Skagit extent skagit historic dates

image

ChristinaB commented 4 years ago

Here is code for processing binary grid outputs into ascii files for Skagit Model Skagit List of Max Annual Dates - continuous 1960-2011 and from 2006-2099 for future data

Number of Rows = 1020 # Number of rows Number of Columns = 916 # Number of columns

Future 95 year for depth to groundwater, 61 years for ice and snow ./myconvert float ascii Map.Soil.TableDepth.bin Map.Soil.TableDepth.asc 96900 916 ./myconvert float ascii Map.Snow.Swq.bin Map.Snow.Swq.asc 62220 916 ./myconvert float ascii Map.Snow.Iwq.bin Map.Snow.Iwq.asc 62220 916

Historic 52 year for depth to groundwater, 43 years for ice and snow ./myconvert float ascii Map.Soil.TableDepth.bin Map.Soil.TableDepth.asc 53040 916 ./myconvert float ascii Map.Snow.Swq.bin Map.Snow.Swq.asc 43860 916 ./myconvert float ascii Map.Snow.Iwq.bin Map.Snow.Iwq.asc 43860 916

ChristinaB commented 4 years ago

Script for converting binary to ascii

foreach F(output_*) cd $F ./myconvert float ascii Map.Soil.TableDepth.bin Map.Soil.TableDepth.asc 96900 916 ./myconvert float ascii Map.Snow.Swq.bin Map.Snow.Swq.asc 62220 916 ./myconvert float ascii Map.Snow.Iwq.bin Map.Snow.Iwq.asc 62220 916 cd .. end

ChristinaB commented 4 years ago
ChristinaB commented 4 years ago

foreach F(output_SCLlandslide_bc*) cp $F/Map.Soil.TableDepth.asc ../dtw_20200320/$F cp $F/Map.Snow.Swq.asc ../ice_20200320/$F cp $F/Map.Snow.Iwq.asc ../snow_20200320/$F end

ChristinaB commented 4 years ago
ChristinaB commented 4 years ago