updates

eviatarbach · Mar 6, 2024 · 3266d6b · 3266d6b
1 parent aa279d8
commit 3266d6b
Show file tree

Hide file tree

Showing 11 changed files with 21,421 additions and 14,067 deletions.
diff --git a/collapse_02.py b/collapse_02.py
@@ -5,7 +5,7 @@
 n_days = 275  # May, June, July, August, September
 D = 4964
 
-prec = xarray.open_dataset('anomalies.nc')
+prec = xarray.open_dataset('data/anomalies.nc')
 
 m = prec["p"]["time"].dt.month
 mask = (m >= 3) & (m <= 11)
@@ -20,4 +20,4 @@
 for i, year in enumerate(years):
     prec_full[:, :, i] = prec_masked[prec_masked.time.dt.year == year].values[:, mask_latlon]
 
-xarray.DataArray(prec_full, name="prec").to_netcdf('prec_out.nc')
+xarray.DataArray(prec_full, name="prec").to_netcdf('data/prec_out.nc')
diff --git a/enoc.ipynb b/enoc.ipynb
@@ -3382,9 +3382,21 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 397,
+   "execution_count": 1,
    "metadata": {},
-   "outputs": [],
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'run_pc_fcst' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "Input \u001b[0;32mIn [1]\u001b[0m, in \u001b[0;36m<cell line: 12>\u001b[0;34m()\u001b[0m\n\u001b[1;32m     15\u001b[0m fcst_pcs \u001b[38;5;241m=\u001b[39m []\n\u001b[1;32m     16\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m year \u001b[38;5;129;01min\u001b[39;00m \u001b[38;5;28mrange\u001b[39m(\u001b[38;5;241m1993\u001b[39m, \u001b[38;5;241m2016\u001b[39m):\n\u001b[0;32m---> 17\u001b[0m     \u001b[43mrun_pc_fcst\u001b[49m(year, month, lead_time, m_p)\n\u001b[1;32m     18\u001b[0m all_true\u001b[38;5;241m.\u001b[39mappend(true_pcs)\n\u001b[1;32m     19\u001b[0m all_fcst\u001b[38;5;241m.\u001b[39mappend(fcst_pcs)\n",
+      "\u001b[0;31mNameError\u001b[0m: name 'run_pc_fcst' is not defined"
+     ]
+    }
+   ],
    "source": [
     "pc_radius = 1.0\n",
     "\n",
@@ -23854,7 +23866,7 @@
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "Python 3",
+   "display_name": "Python 3 (ipykernel)",
    "language": "python",
    "name": "python3"
   },
@@ -23868,7 +23880,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.8.5"
+   "version": "3.9.15"
   }
  },
  "nbformat": 4,

diff --git a/enoc_new.ipynb b/enoc_new.ipynb
diff --git a/mean_anomalies_01.py b/mean_anomalies_01.py
@@ -1,12 +1,12 @@
 import xarray
 
-obs = xarray.open_dataset('prec.nc')
+obs = xarray.open_dataset('data/prec_2021.nc')
 clim = obs["p"].loc[(obs.time.dt.month <= 11) & (obs.time.dt.month >= 3) & (obs.time.dt.year <= 2007) & (obs.time.dt.year >= 1950)]
 clim_grouped = clim.assign_coords(dayofyear=clim.time.dt.strftime("%m-%d")).groupby('dayofyear')
 clim_mean = clim_grouped.mean()
-clim_mean.to_netcdf("clim_mean.nc")
+clim_mean.to_netcdf("data/clim_mean.nc")
 
-all_dat = obs["p"].loc[(obs.time.dt.month <= 11) & (obs.time.dt.month >= 3) & (obs.time.dt.year <= 2007)]
+all_dat = obs["p"].loc[(obs.time.dt.month <= 11) & (obs.time.dt.month >= 3)]
 all_grouped = all_dat.assign_coords(dayofyear=all_dat.time.dt.strftime("%m-%d")).groupby('dayofyear')
 anomalies = all_grouped - clim_mean
-anomalies.to_netcdf('anomalies.nc')
+anomalies.to_netcdf('data/anomalies.nc')
diff --git a/miso_eof_04.jl b/miso_eof_04.jl
@@ -4,13 +4,14 @@ using Statistics
 using HDF5
 
 D = 4964
-N = 31625
 n_pcs = 2
 
 r = h5read("r_summed.h5", "r_summed")
 r = permutedims(r, [2, 1, 3])
 r = reshape(r, D, :)'
 
+N = size(r)[1]
+
 C = r'*r/N
 u, s, v = svd(C)
 pcs = r*v[:, 1:n_pcs]
@@ -20,3 +21,10 @@ eofs = v[:, 1:n_pcs] ./ std_pcs
 
 h5write("eofs.h5", "eofs", eofs)
 h5write("pcs.h5", "pcs", pcs)
+
+r = h5read("r_summed_fcst.h5", "r_summed")
+r = permutedims(r, [2, 1, 3])
+r = reshape(r, D, :)'
+
+pcs = r*eofs
+h5write("pcs_fcst.h5", "pcs", pcs)
diff --git a/prepare_nn_era_centered.py b/prepare_nn_era_centered.py
@@ -10,11 +10,11 @@
 D = 1190
 year0 = 1950
 
-anomalies = xarray.open_dataarray('era5_anomalies.nc')
+anomalies = xarray.open_dataarray('data/era5_anomalies.nc')
 x = anomalies.sel(time=anomalies.time.dt.year <= 1992).stack(stacked=['lat', 'lon']).values
 x = x.reshape(n_years_total, n_days, D)
 pcs = xarray.open_dataset('pcs.h5')["pcs"].values.T
-pcs = pcs.reshape(115, n_days, n_pcs)
+pcs = pcs.reshape(107, n_days, n_pcs)
 
 n_samples = n_years_training*(n_days - 122)
 

diff --git a/prepare_nn_era_future.py b/prepare_nn_era_future.py
@@ -10,11 +10,11 @@
 D = 1190
 year0 = 1950
 
-anomalies = xarray.open_dataarray('era5_anomalies.nc')
+anomalies = xarray.open_dataarray('data/era5_anomalies.nc')
 x = anomalies.sel(time=anomalies.time.dt.year <= 1992).stack(stacked=['lat', 'lon']).values
 x = x.reshape(n_years_total, n_days, D)
 pcs = xarray.open_dataset('pcs.h5')["pcs"].values.T
-pcs = pcs.reshape(115, n_days, n_pcs)
+pcs = pcs.reshape(107, n_days, n_pcs)
 
 n_samples = n_years_training*(n_days - 122)
 

diff --git a/process_era.ipynb b/process_era.ipynb
@@ -792,7 +792,7 @@
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "Python 3",
+   "display_name": "Python 3 (ipykernel)",
    "language": "python",
    "name": "python3"
   },
@@ -806,7 +806,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.8.5"
+   "version": "3.9.15"
   }
  },
  "nbformat": 4,

diff --git a/regress_nn_centered_era.py b/regress_nn_centered_era.py
@@ -29,10 +29,10 @@
 model = Model(inputs=input_layer, outputs=output)
 model.compile(loss="mean_squared_error", optimizer="adam", metrics=["mean_squared_error"])
 
-history = model.fit(X_train, y_train, epochs=200, verbose=1)
+history = model.fit(X_train, y_train, epochs=500, verbose=1)
 
-pickle.dump(sc, open("scaler_centered_era", "wb"))
+pickle.dump(sc, open("scaler_centered_era2", "wb"))
 
-model.save('nn_centered_era')
+model.save('nn_centered_era2')
 
 #print(model.evaluate(X_test, y_test))
diff --git a/regress_nn_future_era.py b/regress_nn_future_era.py
@@ -29,7 +29,7 @@
 model = Model(inputs=input_layer, outputs=output)
 model.compile(loss="mean_squared_error", optimizer="adam", metrics=["mean_squared_error"])
 
-history = model.fit(X_train, y_train, epochs=200, verbose=1)
+history = model.fit(X_train, y_train, epochs=500, verbose=1)
 
 pickle.dump(sc, open("scaler_future_era", "wb"))
 

diff --git a/ssa_prec_03.jl b/ssa_prec_03.jl
@@ -4,15 +4,20 @@ using Statistics
 using NetCDF
 using HDF5
 
-include("../ssa.jl")
+include("ssa.jl")
 using .SSA
 
 M = 61
+n_years = 107
 
-x = ncread("prec_out.nc", "prec")
+x = ncread("data/prec_out.nc", "prec")
 x = permutedims(x, [3, 2, 1])
+x = Array{Float32, 3}(x)
+
+x_ssa = x[:, :, 1:n_years]
+
+ssa_info = SSA.ssa_decompose(x_ssa, M)
 
-ssa_info = SSA.ssa_decompose(x, M)
 r = SSA.ssa_reconstruct(ssa_info, 1:2)
 r_summed = sum(r, dims=1)[1, :, :, :]
 
@@ -22,3 +27,15 @@ h5write("X.h5", "X", ssa_info.X)
 h5write("C.h5", "C", ssa_info.C)
 h5write("r.h5", "r", r)
 h5write("r_summed.h5", "r_summed", r_summed)
+
+x_ssa_fcst = x[:, :, n_years+1:end]
+ssa_info = SSA.ssa_decompose(x_ssa_fcst, M)
+
+ssa_info.eig_vals = eig_vals
+ssa_info.eig_vecs = eig_vecs
+
+r = SSA.ssa_reconstruct(ssa_info, 1:2)
+r_summed = sum(r, dims=1)[1, :, :, :]
+
+h5write("r_fcst.h5", "r", r)
+h5write("r_summed_fcst.h5", "r_summed", r_summed)