<jupyter_start><jupyter_text>Marine Heatwaves in OISST Datavia  [`marineHeatWaves.py`](https://github.com/ecjoliver/marineHeatWaves/blob/master/marineHeatWaves.py)<jupyter_code>fs = gcsfs.GCSFileSystem(project=GCP_PROJECT_ID, token="../gc-<PASSWORD>")
oisst = xr.open_zarr(fs.get_mapper(OISST_GCP))
oisst = oisst.assign_coords(lon=(((oisst.lon + 180) % 360) - 180)).sortby('lon')
oisst
oisst.sst.data
oisst
# s
PNW_LAT = slice(40, 50)
PNW_LON = slice(-165.9, -121.9)
oisst_pnw = oisst.sel(lat = PNW_LAT, lon = PNW_LON).compute()
oisst_pnw.sst.data

oisst_pnw.sst.mean('time').plot()
oisst_pnw.sst.mean(['lat', 'lon']).plot()
def mhw_1d(temps, time):
    if(np.isnan(temps).any()): return xr.DataArray(
        np.full_like(temps, np.nan), 
        name = 'mhw', 
        dims = 'time', 
        coords ={
            'time': time
        }
    )

    ordinals =  np.array([Timestamp(t).toordinal() for t in time])
    dets = mh.detect(ordinals, temps.copy())
    events = dets[0]['n_events']
    del dets[0]['n_events']
    
    arrays = {
        key: xr.DataArray(np.full_like(time, 0., dtype='float64'),
                          name = key, dims='time',
                          coords={
                              'time' : time,
                          })
        for key in ['mhw_cum_intensity']
    }
    
    for event in range(events):
        start_date = dets[0]['date_start'][event]
        end_date = dets[0]['date_end'][event]
        cum_intense = dets[0]['intensity_cumulative'][event]

#         arrays['mhw'].loc[dict(time=slice(start_date, end_date))] = True
        arrays['mhw_cum_intensity'].loc[dict(time=slice(start_date, end_date))] = cum_intense
        
#     clim = xr.DataArray(dets[1]['thresh'], name='clim_thresh', dims='time', coords={'time': time})
        
    return (arrays['mhw_cum_intensity'])
a = xr.apply_ufunc(
    mhw_1d, 
    oisst_pnw.sst.chunk({'lat': 1, 'lon': 1, 'time': -1}), 
    oisst_pnw.time, 
    input_core_dims = [['time'], ['time']],
    output_core_dims=[["time"]], 
    vectorize=True, 
    dask='parallelized', 
    output_dtypes=['float64']
)
a.data

ans = a.compute()
ans.to_netcdf("pnw_mhw_intensity.nc")
sels = ans.isel(lat=25, lon=20)
fig = plt.figure(figsize=(15, 5))

plt.vlines(sels.time[sels].values, 0, 10, color='red', )
plt.plot()
sels = ans.isel(lat=5, lon=20)

import matplotlib.transforms as mt

fig, ax = plt.subplots(figsize=(15, 5))
sst = oisst_pnw.sel(lat=sels.lat.values, lon=sels.lon.values).persist().sst 
ax.vlines(sels.time[sels].values, 0, 10, color='red', alpha=0.2,  transform = trans)

ax.plot(sst.time, sst)
trans = mt.blended_transform_factory(ax.transData, ax.transAxes)





sst = oisst_pnw.sel(lat=sels.lat.values, lon=sels.lon.values).sst 
plt.plot(sst.time, sst)
se
sels
sel
# s
PNW_LAT = slice(44.5, 46.5)
PNW_LON = slice(165.5, 166)
oisst_pnw = oisst.sel(lat = PNW_LAT, lon = PNW_LON).load()
oisst_pnw
a = xr.apply_ufunc(
    mhw_1d, 
    oisst_pnw.sst.chunk({'lat': 5, 'lon': 5, 'time': -1}).load(), 
    oisst_pnw.time.load(), 
    input_core_dims = [['time'], ['time']],
    output_core_dims=[["time"], ["time"]], 
    vectorize=True, 
    output_dtypes=['bool']
)
ans = a.compute()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1812247435498912, 'fraction_numerical': 0.023313646254274168, 'mean_word_length': 3.3310901749663526, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 6, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '1', 'dir': 'jupyter-structured-clean-dedup', 'id': '325275', 'license': 'MIT', 'n_tokens_mistral': 1409, 'n_tokens_neox': 1328, 'n_words': 199, 'path': 'mhw_pipeline/working-oisst-mhw-example.ipynb', 'repo_name': 'HuckleyLab/phyto-mhw', 'size': '3218'}