Skip to content

Commit 41f4a5f

Browse files
committed
updated exposure func. Calls the tidal modelling routine a single time only which streamlines the function
1 parent d0ec967 commit 41f4a5f

File tree

2 files changed

+533
-1198
lines changed

2 files changed

+533
-1198
lines changed

intertidal/exposure.py

Lines changed: 67 additions & 163 deletions
Original file line numberDiff line numberDiff line change
@@ -102,7 +102,7 @@ def temporal_filters(x,
102102
tz=pytz.UTC)
103103

104104
# Replace the UTC datetimes from all_timerange with local times
105-
modelledtides = pd.DataFrame(index=localtides)
105+
ModTides = pd.DataFrame(index=localtides)
106106

107107
# Return the difference in years for the time-period.
108108
# Round up to ensure all modelledtide datetimes are captured in the solar model
@@ -137,10 +137,10 @@ def temporal_filters(x,
137137

138138
# Remove local timezone timestamp column in modelledtides dataframe. Xarray doesn't handle
139139
# timezone aware datetimeindexes 'from_dataframe' very well.
140-
modelledtides.index = modelledtides.index.tz_localize(tz=None)
140+
ModTides.index = ModTides.index.tz_localize(tz=None)
141141

142142
# Create an xr Dataset from the modelledtides pd.dataframe
143-
mt = modelledtides.to_xarray()
143+
mt = ModTides.to_xarray()
144144

145145
# Filter the modelledtides (mt) by the daytime, nighttime datetimes from the sunriset module
146146
# Modelled tides are designated as either day or night by propogation of the last valid index
@@ -169,9 +169,8 @@ def spatial_filters(
169169
ModelledTides,
170170
timeranges,
171171
calculate_quantiles,
172-
tide_cq_dict,
172+
modelledtides_dict,
173173
dem,
174-
exposure
175174
):
176175

177176
"""
@@ -220,15 +219,15 @@ def spatial_filters(
220219
neappeaks = tide_maxima.isel(time=neap_peaks)
221220
timeranges[str(x)]=pd.to_datetime(neappeaks.time)
222221
# Extract the peak height dates
223-
tide_cq = neappeaks.quantile(q=calculate_quantiles,dim='time')
222+
modelledtides = neappeaks.quantile(q=calculate_quantiles,dim='time')
224223

225224
if x in ['Spring_high', 'Spring_low']:
226225
# select for indices associated with peaks
227226
springpeaks = modelledtides_flat.isel(time=modelledtides_flat_peaks).to_dataset()
228227
# Save datetimes for calculation of combined filter exposure
229228
timeranges[str(x)]=pd.to_datetime(springpeaks.time)
230229
# Extract the peak height dates
231-
tide_cq = springpeaks.quantile(q=calculate_quantiles,dim='time')
230+
modelledtides = springpeaks.quantile(q=calculate_quantiles,dim='time')
232231

233232
if x == 'Hightide':
234233
# calculate all the high tide maxima
@@ -244,7 +243,7 @@ def spatial_filters(
244243
# should have approximately equal proportions of daytime and nighttime hightide peaks
245244
if len(lowhigh_peaks)/len(high_peaks) < 0.2:
246245
timeranges[str(x)] = pd.to_datetime(high_peaks2.time)
247-
tide_cq = high_peaks2.quantile(q=calculate_quantiles,dim='time').to_dataset()
246+
modelledtides = high_peaks2.quantile(q=calculate_quantiles,dim='time').to_dataset()
248247
else:
249248
# interpolate the lower hightide curve
250249
low_high_linear = interp(np.arange(0,len(modelledtides_flat)),
@@ -254,7 +253,7 @@ def spatial_filters(
254253
hightide = modelledtides_flat.where(modelledtides_flat >= low_high_linear, drop=True)
255254
## Save datetimes for calculation of combined filter exposure
256255
timeranges[str(x)] = pd.to_datetime(hightide.time)
257-
tide_cq = hightide.quantile(q=calculate_quantiles,dim='time').to_dataset()
256+
modelledtides = hightide.quantile(q=calculate_quantiles,dim='time').to_dataset()
258257

259258
if x == 'Lowtide':
260259
# calculate all the low tide maxima
@@ -270,7 +269,7 @@ def spatial_filters(
270269
# should have approximately equal proportions of daytime and nighttime lowtide peaks
271270
if len(highlow_peaks)/len(low_peaks) < 0.2:
272271
timeranges[str(x)] = pd.to_datetime(low_peaks2.time)
273-
tide_cq = low_peaks2.quantile(q=calculate_quantiles,dim='time').to_dataset()
272+
modelledtides = low_peaks2.quantile(q=calculate_quantiles,dim='time').to_dataset()
274273
else:
275274
# interpolate the higher lowtide curve
276275
high_low_linear = interp(np.arange(0,len(modelledtides_flat)),
@@ -280,20 +279,12 @@ def spatial_filters(
280279
lowtide = modelledtides_flat.where(modelledtides_flat <= high_low_linear, drop=True)
281280
## Save datetimes for calculation of combined filter exposure
282281
timeranges[str(x)] = pd.to_datetime(lowtide.time)
283-
tide_cq = lowtide.quantile(q=calculate_quantiles,dim='time').to_dataset()
284-
285-
# Add tide_cq to output dict
286-
tide_cq_dict[str(x)]=tide_cq.tide_m
287-
# Calculate the tide-height difference between the elevation value and
288-
# each percentile value per pixel
289-
diff = abs(tide_cq.tide_m - dem)
290-
# Take the percentile of the smallest tide-height difference as the
291-
# exposure % per pixel
292-
idxmin = diff.idxmin(dim="quantile")
293-
# Convert to percentage
294-
exposure[str(x)] = idxmin * 100
282+
modelledtides = lowtide.quantile(q=calculate_quantiles,dim='time').to_dataset()
283+
284+
# Add modelledtides to output dict
285+
modelledtides_dict[str(x)]=modelledtides.tide_m
295286

296-
return timeranges, tide_cq_dict, exposure
287+
return timeranges, modelledtides_dict#, exposure
297288

298289
def exposure(
299290
dem,
@@ -317,9 +308,16 @@ def exposure(
317308
the elevation value and modelled tide height percentiles.
318309
319310
For an 'unfiltered', all of epoch-time, analysis, exposure is
320-
calculated per pixel. All of the filter options calculate
311+
calculated per pixel. All other filter options calculate
321312
exposure from a high temporal resolution tide model that is generated
322313
for the center of the nominated area of interest only.
314+
315+
This function firstly calculates a high temporal resolution
316+
tidal model for area (or pixels) of interest. Filtered datetimes and
317+
associated tide heights are then isolated from the tidal model.
318+
Exposure is calculated by comparing the quantiled distribution curve
319+
of modelled tide heights from the filtered datetime dataset with dem
320+
pixel elevations to identify exposure %.
323321
324322
Parameters
325323
----------
@@ -402,7 +400,7 @@ def exposure(
402400
A dictionary of xarray.Datasets containing a named exposure dataset for each
403401
nominated filter, representing the percentage time exposurs of each pixel from seawater
404402
for the duration of the associated filtered time period between `start` and `end`.
405-
tide_cq : dict
403+
modelledtides : dict
406404
A dictionary of xarray.Datasets containing a named dataset of the quantiled high temporal
407405
resolution tide modelling for each filter. Dimesions should be
408406
'quantile', 'x' and 'y'.
@@ -449,7 +447,7 @@ def exposure(
449447
## Create empty datasets to store outputs into
450448
exposure = xr.Dataset(coords=dict(y=(['y'], dem.y.values),
451449
x=(['x'], dem.x.values)))
452-
tide_cq_dict = xr.Dataset(coords=dict(y=(['y'], dem.y.values),
450+
modelledtides_dict = xr.Dataset(coords=dict(y=(['y'], dem.y.values),
453451
x=(['x'], dem.x.values)))
454452

455453
## Create an empty dict to store temporal `time_range` variables into
@@ -472,72 +470,31 @@ def exposure(
472470

473471
# Calculate a tidal model. Run at pixel resolution if any filter is 'unfiltered' else run at low res
474472
if 'unfiltered' in filters:
475-
tide_cq, _ = pixel_tides_ensemble(
473+
mod_tides, _ = pixel_tides_ensemble(
476474
dem,
477475
model=tide_model,
478476
calculate_quantiles=calculate_quantiles,
479477
times=time_range,
480478
directory=tide_model_dir,
481479
ancillary_points="data/raw/tide_correlations_2017-2019.geojson",
482480
)
483-
tide_cq_flat = tide_cq.mean(dim=["x","y"])
484-
else:
485-
#Calculate a low spatial res tidal model
486-
tide_cq = pixel_tides_ensemble(
487-
dem,
488-
model=tide_model,
489-
# calculate_quantiles=calculate_quantiles,
490-
times=time_range,
491-
directory=tide_model_dir,
492-
ancillary_points="data/raw/tide_correlations_2017-2019.geojson",
493-
resample=False
494-
)
495-
496-
# Flatten low res tidal model
497-
## To reduce compute, average across the y and x dimensions
498-
tide_cq_flat = tide_cq.mean(dim=["x","y"])
499-
500-
# Calculate exposure using pixel-based tide modelling for unfiltered, all of epoch time period
501-
if 'unfiltered' in filters:
502-
print('-----\nCalculating unfiltered exposure')
503-
504-
# tide_cq, _ = pixel_tides_ensemble(
505-
# dem,
506-
# model=tide_model,
507-
# calculate_quantiles=calculate_quantiles,
508-
# times=time_range,
509-
# directory=tide_model_dir,
510-
# ancillary_points="data/raw/tide_correlations_2017-2019.geojson",
511-
# )
512-
513-
# Add tide_cq to output dict
514-
tide_cq_dict['unfiltered']=tide_cq
515-
# Calculate the tide-height difference between the elevation value and
516-
# each percentile value per pixel
517-
diff = abs(tide_cq - dem)
518-
# Take the percentile of the smallest tide-height difference as the
519-
# exposure % per pixel
520-
idxmin = diff.idxmin(dim="quantile")
521-
# Convert to percentage
522-
exposure['unfiltered'] = idxmin * 100
523-
524-
# Prepare for spatial filtering. Calculate the pixel-based all-epoch high res tide model.
525-
# Reduce the tide-model to the mean for the area of interest (reduce compute).
526-
# if any (x in sptl_filters for x in filters):
527-
# print ('-----\nCalculating tide model for spatial filters')
528-
# # Use low res modelling for spatial filters
529-
# ModelledTides = pixel_tides_ensemble(
530-
# dem,
531-
# times=time_range,
532-
# directory=tide_model_dir,
533-
# model = tide_model,
534-
# ancillary_points="data/raw/tide_correlations_2017-2019.geojson",
535-
# resample=False
536-
# )
537-
538-
# ## To reduce compute, average across the y and x dimensions
539-
# modelledtides_flat = ModelledTides.mean(dim=["x","y"])
540-
481+
# Add modelledtides to output dict
482+
modelledtides_dict['unfiltered']=mod_tides
483+
484+
# For all other filter types, calculate a low spatial res tidal model
485+
modelledtides = pixel_tides_ensemble(
486+
dem,
487+
model=tide_model,
488+
times=time_range,
489+
directory=tide_model_dir,
490+
ancillary_points="data/raw/tide_correlations_2017-2019.geojson",
491+
resample=False
492+
)
493+
494+
# Flatten low res tidal model
495+
## To reduce compute, average across the y and x dimensions
496+
modelledtides_flat = modelledtides.mean(dim=["x","y"])
497+
541498
# Filter the input timerange to include only dates or tide ranges of interest
542499
# if filters is not None:
543500
for x in filters:
@@ -550,103 +507,50 @@ def exposure(
550507
dem)
551508

552509
elif x in sptl_filters:
553-
print(f'-----\nCalculating {x} exposure')
510+
print(f'-----\nCalculating {x} timerange')
554511

555-
timeranges, tide_cq_dict, exposure = spatial_filters(
556-
modelled_freq,
557-
x,
558-
tide_cq_flat,
559-
tide_cq,
560-
timeranges,
561-
calculate_quantiles,
562-
tide_cq_dict,
563-
dem,
564-
exposure
565-
)
512+
timeranges, modelledtides_dict = spatial_filters(
513+
modelled_freq,
514+
x,
515+
modelledtides_flat,
516+
modelledtides,
517+
timeranges,
518+
calculate_quantiles,
519+
modelledtides_dict,
520+
dem,
521+
)
566522

567-
## Intersect the filters of interest to extract the common datetimes for calc of combined filters
523+
## Intersect the filters of interest to extract the common datetimes for calculation of combined filters
568524
if filters_combined is not None:
569525
for x in filters_combined:
570526
y=x[0]
571527
z=x[1]
572528
timeranges[str(y+"_"+z)] = timeranges[y].intersection(timeranges[z])
573529

574-
## Generator expression to calculate exposure for each nominated filter in temp_filters
530+
# Intersect datetimes of interest with the low-res tidal model
575531
# Don't calculate exposure for spatial filters. This has already been calculated.
576-
gen = (x for x in timeranges if x not in sptl_filters)
577-
'''
578-
TEsting new idea
579-
'''
580-
# #Generate single tidal model
581-
# tide_cq = pixel_tides_ensemble(
582-
# dem,
583-
# model=tide_model,
584-
# # calculate_quantiles=calculate_quantiles,
585-
# times=time_range,
586-
# directory=tide_model_dir,
587-
# ancillary_points="data/raw/tide_correlations_2017-2019.geojson",
588-
# resample=False
589-
# )
590-
591-
# # Flatten low res tidal model
592-
# ## To reduce compute, average across the y and x dimensions
593-
# tide_cq_flat = tide_cq.mean(dim=["x","y"])
594-
595-
532+
gen = (x for x in timeranges if x not in sptl_filters)
596533
for x in gen:
597-
print(f'-----\nCalculating {x} exposure')
598534
# Extract filtered datetimes from the full tidal model
599-
tide_cq_flat = tide_cq_flat.sel(time=timeranges[str(x)])
535+
modelledtides_x = modelledtides_flat.sel(time=timeranges[str(x)])
600536
# Calculate quantile values on remaining tide heights
601-
tide_cq_flat = tide_cq_flat.quantile(q=calculate_quantiles,dim='time').to_dataset().tide_m
602-
# Add tide_cq to output dict
603-
tide_cq_dict[str(x)]=tide_cq_flat
537+
modelledtides_x = modelledtides_x.quantile(q=calculate_quantiles,dim='time').to_dataset().tide_m
538+
# Add modelledtides to output dict
539+
modelledtides_dict[str(x)]=modelledtides_x
540+
541+
# Calculate exposure per filter
542+
for x in modelledtides_dict:
543+
print (f'-----\nCalculating {x} exposure')
604544
# Calculate the tide-height difference between the elevation value and
605545
# each percentile value per pixel
606-
diff = abs(tide_cq_flat - dem)
546+
diff = abs(modelledtides_dict[str(x)] - dem)
607547
# Take the percentile of the smallest tide-height difference as the
608548
# exposure % per pixel
609549
idxmin = diff.idxmin(dim="quantile")
610550
# Convert to percentage
611551
exposure[str(x)] = idxmin * 100
612-
'''
613-
--- end
614-
'''
615-
616-
# for x in gen:
617-
# # Run the pixel_tides function with the calculate_quantiles option and default high-res outputs.
618-
# # For each pixel, an array of tideheights is returned, corresponding
619-
# # to the percentiles from `calculate_quantiles` of the timerange-tide model that
620-
# # each tideheight appears in the model.
621-
622-
# # Print
623-
# print(f'-----\nCalculating {x} exposure')
624-
625-
# tide_cq = pixel_tides_ensemble(
626-
# dem,
627-
# calculate_quantiles=calculate_quantiles,
628-
# times=timeranges[str(x)],
629-
# directory=tide_model_dir,
630-
# ancillary_points="data/raw/tide_correlations_2017-2019.geojson",
631-
# model=tide_model,
632-
# resample=False
633-
# )
634-
# # Flatten low res tidal model
635-
# ## To reduce compute, average across the y and x dimensions
636-
# tide_cq_flat = tide_cq.mean(dim=["x","y"])
637-
638-
# # Add tide_cq to output dict
639-
# tide_cq_dict[str(x)]=tide_cq_flat
640-
# # Calculate the tide-height difference between the elevation value and
641-
# # each percentile value per pixel
642-
# diff = abs(tide_cq_flat - dem)
643-
# # Take the percentile of the smallest tide-height difference as the
644-
# # exposure % per pixel
645-
# idxmin = diff.idxmin(dim="quantile")
646-
# # Convert to percentage
647-
# exposure[str(x)] = idxmin * 100
648-
649-
return exposure, tide_cq_dict
552+
553+
return exposure, modelledtides_dict
650554

651555

652556

0 commit comments

Comments
 (0)