run_swr.py 2.89 KB
Newer Older
1
import os
2
3
4
5
import pandas as pd
import numpy as np
import glob
import re
6
7
import datetime

8
9
10
11
12
from scipy.interpolate import interp1d

from utils.sunposition import sunpos
import utils.utils as u
import utils.auxdata as ua
13
from trios.process import *
14
15
16
17
18
19
20

plot=False
odir = os.path.abspath('/DATA/OBS2CO/data/trios/surface_water')
dirfig = os.path.abspath('/DATA/OBS2CO/data/trios/fig')

swrfiles = glob.glob("/DATA/OBS2CO/data/trios/raw/Lu0*idpr*.csv")

21
22
23
24
25
26
27
28
coordf = glob.glob("/DATA/OBS2CO/data/info/mesures_in_situ_test.csv")[0]
coords = pd.read_csv(coordf, sep=';')
coords['Date_prel']=pd.to_datetime(coords['Date_prel'])
coords['h_debut']=coords['Date_prel'] + pd.to_timedelta(coords['h_debut'])
coords['h_fin']=coords['Date_prel'] + pd.to_timedelta(coords['h_fin'])
# TODO warning: time is set as start_time + 15 minutes (can be more accurate)
coords['Date_prel']= coords['h_debut']+datetime.timedelta(minutes = 15)

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
iopw = ua.iopw()
iopw.load_iopw()

def add_curve(ax, x, mean, std, c='red', label=''):
    ax.plot(x, mean, linestyle='solid', c=c, lw=2.5,
            alpha=0.8, label=label)
    ax.fill_between(x,
                    mean - std,
                    mean + std, alpha=0.35, color=c)


idpr = '167'

# get idpr numbers
idprs = np.unique([re.findall(r'idpr(\d+)', x)[0] for x in swrfiles])
#idprs = np.array(['170'])
# loop over idpr
for idpr in idprs:
    c = coords[coords.ID_prel == int(idpr)]  # .values[0]
    date=c['Date_prel'].dt.strftime('%Y%m%d')
    lat = c['Lat'].values[0]
    lon = c['Lon'].values[0]
    alt = 0  # c['Altitude']
    name = c['ID_lac'].values[0]

54
55
    ofile=os.path.join(odir,'Rrs_swr_'+date.values[0]+'_idpr'+idpr+'_'+name+'.csv')
    header = c.stack(dropna=False)
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
    header.index = header.index.droplevel()
    header.to_csv(ofile, header=None)

    # -----------------------------------------------
    #   SWR processing
    # -----------------------------------------------

    uswr = u.swr_data(idpr, swrfiles)
    if uswr.file:
        df, wl_swr = uswr.reader(lat, lon, alt)
        df['sza', ''] = np.nan
        for index, row in df.iterrows():
            # print index
            sza = sunpos(index, lat, lon, alt)[1]
            df.at[index, 'sza'] = sza
        swr = swr_process(df, wl_swr)
        Rrs_swr = swr.call_process()
        Rrs_stat = Rrs_swr.describe()
        Rrs_stat.columns=Rrs_stat.columns.droplevel()
        Rrs_stat = Rrs_stat.T
        Rrs_stat.to_csv(ofile,mode='a')
        if plot:
78
79
            import matplotlib.pyplot as plt

80
81
82
83
84
            fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(10, 8))
            fig.subplots_adjust(left=0.1, right=0.9, hspace=.5, wspace=0.65)
            add_curve(ax, wl_swr, Rrs_swr.transpose().mean(axis=1), Rrs_swr.transpose().std(axis=1), label='swr', c='black')
            Rrs_swr = swr.call_process(shade_corr=True)
            add_curve(ax, wl_swr, Rrs_swr.transpose().mean(axis=1), Rrs_swr.transpose().std(axis=1), label='swr', c='red')