forked from lgloege/bacalhau_socat_test
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathmain.py
78 lines (59 loc) · 2.3 KB
/
main.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
import fsspec
import xarray as xr
import pandas as pd
import numpy as np
import pyseaflux
def lon_360_to_180(ds=None, lonVar=None):
lonVar = "lon" if lonVar is None else lonVar
return (ds.assign_coords({lonVar: (((ds[lonVar] + 180) % 360) - 180)})
.sortby(lonVar)
.astype(dtype='float32', order='C'))
def center_dates(ds):
# start and end date
start_date = str(ds.time[0].dt.strftime('%Y-%m').values)
end_date = str(ds.time[-1].dt.strftime('%Y-%m').values)
# monthly dates centered on 15th of each month
dates = pd.date_range(start=f'{start_date}-01T00:00:00.000000000',
end=f'{end_date}-01T00:00:00.000000000',
freq='MS') + np.timedelta64(14, 'D')
return ds.assign(time=dates)
def get_and_process_sst(url=None):
# get noaa sst
if url is None:
url = ("input/sst.mnmean.nc")
with fsspec.open(url) as fp:
ds = xr.open_dataset(fp)
ds = lon_360_to_180(ds)
ds = center_dates(ds)
return ds
def get_and_process_socat(url=None):
if url is None:
url = ("input/SOCATv2021_tracks_gridded_monthly.nc.zip")
with fsspec.open(url, compression='zip') as fp:
ds = xr.open_dataset(fp)
ds = ds.rename({"xlon": "lon", "ylat": "lat", "tmnth": "time"})
ds = center_dates(ds)
return ds
def main():
# get sst and socat
ds_sst = get_and_process_sst()
ds_socat = get_and_process_socat()
# merge datasets together
time_slice = slice("1981-12", "2022-05")
ds_out = xr.merge([ds_sst['sst'].sel(time=time_slice),
ds_socat['fco2_ave_unwtd'].sel(time=time_slice)])
# calculate pco2 from fco2
ds_out['pco2_ave_unwtd'] = xr.apply_ufunc(
pyseaflux.fCO2_to_pCO2,
ds_out['fco2_ave_unwtd'],
ds_out['sst'])
# include metadata
ds_out['pco2_ave_unwtd'].attrs['units'] = 'uatm'
ds_out['pco2_ave_unwtd'].attrs['notes'] = ("calculated using" +
"NOAA OI SST V2" +
"and pyseaflux package")
# save as zarr
ds_out.to_zarr("output/SOCATv2021_tracks_gridded_monthly_processed.zarr")
print("Zarr file written to disk, job completed successfully")
if __name__ == "__main__":
main()