From 8a261f81749d17772e0ba49d78bea53462592769 Mon Sep 17 00:00:00 2001 From: Muaz Ahmad Date: Wed, 2 Aug 2023 12:08:51 +0500 Subject: [PATCH] Fixed sst script, leap year handling, proper range, sigma calc --- python_scripts/arc_sst_anom.py | 31 +++++++++++++++++-------------- 1 file changed, 17 insertions(+), 14 deletions(-) diff --git a/python_scripts/arc_sst_anom.py b/python_scripts/arc_sst_anom.py index 8a8b40e..f62fbf2 100644 --- a/python_scripts/arc_sst_anom.py +++ b/python_scripts/arc_sst_anom.py @@ -3,22 +3,25 @@ import matplotlib.pyplot as plt import numpy, pandas import requests -res = requests.get('https://climatereanalyzer.org/clim/sst_daily/json/oisst2.1_world2_sst_day.json') -if res.status_code == 200: - data = json.loads(res.text) +start_range = 29 +start_year = 1982 +res = requests.get('https://climatereanalyzer.org/clim/sst_daily/json/oisst2.1_natlan1_sst_day.json') +data = json.loads(res.text) +data = data[1:-3] -data = data[:-3] # remove final three cols (avg and sds) - -temps = numpy.array([i['data'] for i in data], numpy.float32).T -cols = [i['name'] for i in data] +temps = numpy.array([i['data'] for i in data], numpy.float32).T +cols = [i['name'] for i in data] df = pandas.DataFrame(temps, columns = cols) -diffs = pandas.DataFrame(columns = cols, dtype=numpy.float32) -for col in cols[20:]: - i = cols.index(col) - means = df[df.columns[0:20]].mean(axis=1) # 1981 - 2001 mean sst - diffs[col] = df[col] - means -series = diffs.melt().drop('variable', axis=1) +mean = df[df.columns[:start_range + 1]].mean(axis=1) +sd = df[df.columns[:start_range + 1]].std(axis=1) +diffs = pandas.DataFrame(columns = cols, dtype=numpy.float32) +for col in cols[start_range:]: + if int(col) % 4 != 0: + df[col].iloc[60:] = df[col].iloc[60:].shift(1) + i = cols.index(col) + diffs[col] = (df[col] - mean) / sd +series = diffs.melt().drop('variable', axis=1) series.plot(figsize=(20,5)) -plt.xticks(numpy.arange(0,len(series),366)[20:], (numpy.arange(0, len(series), 366)[20:] / 366 + 1981).astype(int)) +plt.xticks(numpy.arange(0,len(series),366)[start_range+1:], (numpy.arange(0, len(series), 366)[start_range+1:] / 366 + start_year).astype(int)) # fix labels plt.show()