run_swr.py 3.11 KB
Newer Older
1
2
import glob
import re
3
4
import datetime

5
6
from trios.utils.sunposition import sunpos
from trios.utils import utils as u
7
from trios.process import *
8

Harmel Tristan's avatar
Harmel Tristan committed
9
plot=True #False
10
odir = os.path.abspath('/DATA/OBS2CO/data/trios/surface_water')
Harmel Tristan's avatar
Harmel Tristan committed
11
dirfig = os.path.join(odir,'fig')
12

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

15
16
17
18
19
20
21
22
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)

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
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]

48
49
    ofile=os.path.join(odir,'Rrs_swr_'+date.values[0]+'_idpr'+idpr+'_'+name+'.csv')
    header = c.stack(dropna=False)
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
    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:
72
            import matplotlib.pyplot as plt
Harmel Tristan's avatar
Harmel Tristan committed
73
            plt.rcParams.update({'font.size': 18})
74
            fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(10, 8))
Harmel Tristan's avatar
Harmel Tristan committed
75
            fig.subplots_adjust(left=0.16, right=0.9, hspace=.5, wspace=0.65)
76
77
78
            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')
Harmel Tristan's avatar
Harmel Tristan committed
79
80
81
82
83
            ax.legend(loc='best', frameon=False)
            ax.set_ylabel(r'$R_{rs}\  (sr^{-1})$')
            ax.set_xlabel(r'$Wavelength (nm)$')
            fig.savefig(os.path.join(dirfig, 'trios_swr_'+date.values[0]+'_'+ name + '_idpr' + idpr + '.png'))
            plt.close()