Merge pull request #17 from OpenGeoscience/new-initial-conditions

New initial conditions - 002 and 003

Author: annehaley

uvdat_flood_sim/__main__.py

@@ -18,6 +18,12 @@ def _ensure_dir_exists(ctx, param, value):
 
 
 @click.command(name='Dynamic Flood Simulation')
+@click.option(
+    '--initial-conditions-id', '-i',
+    type=click.Choice(['001', '002', '003']),
+    default='001',
+    help='Initialization for real-world emissions scenarios to represent future climate conditions',
+)
 @click.option(
     '--time-period', '-t',
     type=click.Choice(['2031-2050', '2041-2060']),
@@ -83,6 +89,7 @@ def _ensure_dir_exists(ctx, param, value):
     help='Library to use for writing result tiff'
 )
 def main(
+    initial_conditions_id: Literal['001', '002', '003'],
     time_period: Literal['2031-2050', '2041-2060'],
     annual_probability: float,
     hydrograph_name: Literal['short_charles', 'long_charles'],
@@ -97,13 +104,14 @@ def main(
     logging.basicConfig(level=logging.INFO)
 
     logger.info((
-        f'Inputs: {time_period=}, {annual_probability=}, {hydrograph=}, '
+        f'Inputs: {initial_conditions_id=}, {time_period=}, {annual_probability=}, {hydrograph=}, '
         f'{pet_percentile=}, {sm_percentile=}, {gw_percentile=}, '
         f'{output_path=}, {animation=}'
     ))
     start = time.perf_counter()
 
     flood = run_sim(
+        initial_conditions_id=initial_conditions_id,
         time_period=time_period,
         annual_probability=annual_probability,
         hydrograph_name=hydrograph_name,
@@ -113,11 +121,18 @@ def main(
         gw_percentile=gw_percentile,
     )
 
-    write_multiframe_geotiff(flood, output_path, writer=tiff_writer)
+    write_multiframe_geotiff(
+        flood_results=flood,
+        output_path=output_path / 'flood_simulation.tif',
+        writer=tiff_writer,
+    )
     logger.info(f'Done in {time.perf_counter() - start} seconds.')
 
     if animation:
-        animate_results(flood, output_path)
+        animate_results(
+            flood_results=flood,
+            output_path=output_path / 'animation.gif',
+        )
 
 
 if __name__ == '__main__':

uvdat_flood_sim/animate_results.py

@@ -4,26 +4,26 @@
 import matplotlib.animation as ani
 
 
-def animate(results, output_folder):
-    n_frames = results.shape[0]
-    vmin, vmax = numpy.min(results), numpy.max(results)
+def animate(*, flood_results, output_path):
+    n_frames = flood_results.shape[0]
+    vmin, vmax = numpy.min(flood_results), numpy.max(flood_results)
 
     fig, ax = plt.subplots()
     im = ax.imshow(
-        results[0],
+        flood_results[0],
         cmap='viridis',
         vmin=vmin,
         vmax=vmax,
     )
 
     def update(i):
         ax.set_title(f'Hour={i + 1}')
-        im.set_data(results[i])
+        im.set_data(flood_results[i])
         return im
 
     animation = ani.FuncAnimation(fig, update, n_frames, interval=1000)
     animation.save(
-        output_folder / 'animation.gif',
+        output_path,
         writer=ani.PillowWriter(fps=2)
     )
 

uvdat_flood_sim/constants.py

@@ -12,8 +12,18 @@
 PERCENTILES_PATH = DATA_FOLDER / 'flood_param_percentiles.json'
 
 CESM_DATA = {
-    '2031-2050': DATA_FOLDER / 'boston_cesm2-001_2031_2050.npy',
-    '2041-2060': DATA_FOLDER / 'boston_cesm2-001_2041_2060.npy',
+    '001': {
+        '2031-2050': DATA_FOLDER / 'boston_cesm2-001_2031_2050.npy',
+        '2041-2060': DATA_FOLDER / 'boston_cesm2-001_2041_2060.npy',
+    },
+    '002': {
+        '2031-2050': DATA_FOLDER / 'boston_cesm2-002_2031_2050.npy',
+        '2041-2060': DATA_FOLDER / 'boston_cesm2-002_2041_2060.npy',
+    },
+    '003': {
+        '2031-2050': DATA_FOLDER / 'boston_cesm2-003_2031_2050.npy',
+        '2041-2060': DATA_FOLDER / 'boston_cesm2-003_2041_2060.npy',
+    },
 }
 
 HYDROGRAPHS = dict(

uvdat_flood_sim/data/boston_cesm2-002_2031_2050.npy

@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:cdcf3d6bdd1353142437d4846b2c353334ced817be2e79c9d928f0b3cb287d99
+size 700928

uvdat_flood_sim/data/boston_cesm2-002_2041_2060.npy

@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:2717b313bbab31da14c7b686a1f9e58999ae06f79e0a533a9e5738c204cf1a47
+size 700928

uvdat_flood_sim/data/boston_cesm2-003_2031_2050.npy

@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:c951f8d0203642a5e73994fe70a4a8e7c42bab8bbd81941e3fec5b2001987ea8
+size 700928

uvdat_flood_sim/data/boston_cesm2-003_2041_2060.npy

@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:3c84df4f798552160b4f3fd6700685d3a8152a8918b04e9942e355f998b1dbf9
+size 700928

uvdat_flood_sim/downscaling_prediction.py

@@ -16,8 +16,9 @@ def annual_precipitation_maxima(daily):
     return apm
 
 
-def downscale_boston_cesm(cesm_id, annual_probability):
-    cesm_file = CESM_DATA[cesm_id] # CESM data is a projection of 20 years of daily weather conditions, spatially coarse over a large region
+def downscale_boston_cesm(initial_conditions_id, time_period, annual_probability):
+    # CESM data is a projection of 20 years of daily weather conditions, spatially coarse over a large region
+    cesm_file = CESM_DATA[initial_conditions_id][time_period]
     cesm_data = numpy.load(cesm_file, allow_pickle=True)
     with DOWNSCALING_MODEL_PATH.open('rb') as m:
         model = pickle.load(m)

uvdat_flood_sim/run.py

@@ -15,6 +15,7 @@
 
 def run_sim(
     *,
+    initial_conditions_id: Literal['001', '002', '003'],
     time_period: Literal['2031-2050', '2041-2060'],
     annual_probability: float,
     hydrograph_name: Literal['short_charles', 'long_charles'] | None = None,
@@ -49,7 +50,7 @@ def run_sim(
     ground_water = percentiles['gw'][gw_percentile]
 
     # Extreme precipitation level in millimeters
-    level = downscale_boston_cesm(time_period, annual_probability)
+    level = downscale_boston_cesm(initial_conditions_id, time_period, annual_probability)
     logger.info(f'Downscaling prediction: precipitation level = {level}')
 
     # Obtain discharge

uvdat_flood_sim/save_results.py

@@ -78,10 +78,9 @@ def large_image_write(results, output_path):
 def write_multiframe_geotiff(
     *,
     flood_results: NDArray[numpy.float32],
-    output_folder: Path,
+    output_path: Path,
     writer: Literal['rasterio', 'large_image'] = 'rasterio',
 ) -> None:
-    output_path = output_folder / 'flood_simulation.tif'
     if writer == 'rasterio':
         rasterio_write(flood_results, output_path)
     elif writer == 'large_image':