En raison du déménagement des baies serveurs, les services gitlab.irstea.fr et mattermost.irstea.fr seront interrompus le samedi 2 octobre 2021 au matin. Ils devraient revenir à la normale dans la journée.

Commit 6b4880fc authored by Harmel Tristan's avatar Harmel Tristan
Browse files

add 'run' file for swr, define output format

parent e46c0a24
import numpy as np
import pandas as pd
import scipy.interpolate as si
from scipy import interpolate, integrate
from scipy.optimize import curve_fit
import plotly.plotly as py
# import plotly.graph_objs as go
from plotly.graph_objs import *
from utils.utils import reshape as r
import utils.auxdata as ua
from config import *
......@@ -30,25 +33,33 @@ class awr_process:
self.rhoM1999.dropna(inplace=True)
self.rhoM2015.dropna(inplace=True)
def get_rho_values(self, sza, vza, azi, ws=2, aot=0.1, wl=[550]):
def get_rho_values(self, sza, vza, azi, ws=[2], aot=[0.1], wl=[550], sunglint=False):
'''
:param sza:
:param vza:
:param azi:
:param ws:
:param aot:
:param wl:
Interpolate the rho factor values from tabulated data
:param sza: solar zenith angle in deg, array-like
:param vza: view zenith angle in deg, array-like
:param azi: relative azimuth in deg (=0 when looking at Sun), array-like
:param ws: wind speed, m/s,(based on Cox-Munk parametrization of surface roughness) array-like
:param aot: aerosol optical thickness at 550 nm, array-like
:param wl: wavelength in nm, array-like
:param sunglint: add sunglint component in rho calculation if True
:return:
'''
grid = self.rho.rho.index.levels
grid = self.rho.rho.index.levels[2:]
# convert pandas dataframe into 6D array of the tabulated rho values for interpolation
rho_ = r().df2ndarray(self.rho, 'rho')
rhoname = 'rho'
if sunglint:
rhoname = 'rho_g'
rho_6d = r().df2ndarray(self.rho, rhoname)
rho_ = calc().spline_2d(grid[-2:], rho_6d, (azi, vza))
rho_wl = calc().spline_4d(grid[:-2], rho_, (ws, aot, wl, sza))
rho_wl = calc().spline_4d(grid, rho_[1, 1, ...], (wl, sza, azi, vza))
return rho_wl.squeeze()
def get_rho_mobley(self, rhodf, sza, vza, azi, ws):
......@@ -78,10 +89,13 @@ class awr_process:
Lt = self.df.loc[:, ("Lt")]
Lsky = self.df.loc[:, ("Lsky")]
Ed = self.df.loc[:, ("Ed")]
sza = self.df.loc[:, ("sza")]
return self.process(wl, Lt, Lsky, Ed, sza, vza, azi, ws, aot)
sza = self.df.loc[:, ("sza")].values.mean()
Rrs = self.process(wl, Lt, Lsky, Ed, sza, vza, azi, ws, aot)
Rrs.columns = pd.MultiIndex.from_product([['Rrs(awr)'], self.Rrs.columns], names=['param', 'wl'])
self.Rrs = Rrs
return self.Rrs
def process(self, wl, Lt, Lsky, Ed, sza, vza=[40], azi=[135], ws=2, aot=0.1):
def process(self, wl, Lt, Lsky, Ed, sza, vza=[40], azi=[135], ws=[2], aot=[0.1]):
'''
:param wl:
......@@ -96,9 +110,8 @@ class awr_process:
:return:
'''
rho = self.get_rho_values([sza.values.mean()],vza,azi,wl=wl )
rho = self.get_rho_values([sza], vza, azi, wl=wl, ws=ws, aot=aot)
self.Rrs = (Lt - rho * Lsky) / Ed
self.Rrs.columns = pd.MultiIndex.from_product([['Rrs(awr)'], self.Rrs.columns], names=['param', 'wl'])
return self.Rrs, rho
......@@ -108,13 +121,78 @@ class swr_process:
self.df = df
self.wl = wl
def process(self):
def call_process(self, shade_corr=False):
wl = self.wl
df = self.df
Rrs = df.loc[:, ("Lu0+")] / df.loc[:, ("Ed")]
Lu = self.df.loc[:, ("Lu0+")]
Ed = self.df.loc[:, ("Ed")]
sza = self.df.loc[:, ("sza")].values.mean()
Rrs = self.process(Lu, Ed, sza, wl, shade_corr=shade_corr)
Rrs.columns = pd.MultiIndex.from_product([['Rrs(swr)'], Rrs.columns], names=['param', 'wl'])
self.Rrs = Rrs
return Rrs
def process(self, Lu, Ed, sza, wl, R=0.05, shade_corr=False):
Rrs = Lu / Ed
ang_w = calc().angle_w(sza)
iopw = ua.iopw()
iopw.load_iopw()
iopw.get_iopw(wl)
a, bb = iopw.aw, iopw.bbw
# TODO add particulate and dissolved component to a and bb values
# a,bb = aux.get_iop(..., withwater=True)
acdom = ua.cdom(0.5, wl).get_acdom()
a = a + acdom + 0.4
bb = bb + 0.05
if shade_corr:
Rrs = self.shade_corr(Rrs, R, ang_w, a, bb, wl)
# Rrs.columns = pd.MultiIndex.from_product([['Rrs(swr)'], Rrs.columns], names=['param', 'wl'])
self.Rrs = Rrs
self.a = a
self.bb = bb
self.acdom = acdom
return self.Rrs
def epsilon(self, K, R, ang_w):
'''
epsilon from Shang et al, 2017, Applied Optics
:param K:
:param R:
:param ang_w: Sun zenith angle below surface (in deg)
:return:
'''
self.eps = np.array(1 - np.exp(-K * R / np.tan(np.radians(ang_w))))
return self.eps
def K(self, a, bb, ang_w):
'''
K (sum attenuation coef. of Lu in and outside the shade) from Shang et al, 2017, Applied Optics
:param a: total absorption coefficient (m-1)
:param bb: total backscattering coefficient (m-1)
:param ang_w: Sun zenith angle below surface (in deg)
:return:
'''
sin_ang_w = np.sin(np.radians(ang_w))
self.K_ = (3.15 * sin_ang_w + 1.15) * a * np.exp(-1.57 * bb) \
+ (5.62 * sin_ang_w - 0.23) * bb * np.exp(-0.5 * a)
return self.K_
def shade_corr(self, Rrs, R, ang_w, a, bb, wl, wl_cutoff=900):
'''
Correction of shading error from Shang et al, 2017, Applied Optics
:param Rrs:
:param R:
:param ang_w:
:param a:
:param bb:
:return:
'''
K = self.K(a, bb, ang_w)
eps = self.epsilon(K, R, ang_w)
eps[wl > wl_cutoff] = 0
self.Rrs = Rrs / (1 - eps)
return self.Rrs
......@@ -127,18 +205,33 @@ class iwr_process:
wl = self.wl
df = self.df
Rrs = df.loc[:, ("Lu0+")] / df.loc[:, ("Ed")]
Rrs.columns = pd.MultiIndex.from_product([['Rrs(swr)'], Rrs.columns], names=['param', 'wl'])
self.Rrs = Rrs
return self.Rrs
reflectance = df.loc[:, ("Luz")] / df.loc[:, ("Edz")]
reflectance.columns = pd.MultiIndex.from_product([['reflectance'], reflectance.columns], names=['param', 'wl'])
self.reflectance = reflectance
df['rounded_depth', ''] = df.prof_Edz.round(1)
df.groupby('rounded_depth').mean()
return self.reflectance
@staticmethod
def f_Edz(depth, Kd, Ed0):
'''simple Edz model for homogeneous water column'''
return Ed0 * np.exp(-Kd*depth)
@staticmethod
def f_logEdz(depth, Kd, Ed0):
'''simple Edz model for homogeneous water column'''
return np.log(1 + iwr_process.f_Edz(depth, Kd, Ed0)) #Ed0) -Kd*depth
def Kd(self, depth, Edz):
Kd = np.diff(Edz) / np.diff(depth)
def plot_raw(self):
def plot_raw(self,x='Luz',y='prof_Luz'):
trace = Scattergl(
x=self.df['Luz'].values,
y=self.df['depth_Luz'].values,
x=self.df[x].values,
y=self.df[y].values,
text=self.df.index.get_level_values(0),
hoverinfo="text",
......@@ -172,8 +265,8 @@ class iwr_process:
title=''
),
yaxis=dict(
range=[min(-5, min(self.df['depth_Luz'])),
max(0, max(self.df['depth_Luz']))],
range=[min(-5, min(self.df[y])),
max(0, max(self.df[y]))],
showline=False,
fixedrange=True,
zeroline=False,
......@@ -200,7 +293,7 @@ class self_shading:
self.eps_dif_EuZ = 2.70
def epsilon(self, sza):
eps = si.interp1d(self.ang, self.eps_dif_EuZ)(sza)
eps = interpolate.interp1d(self.ang, self.eps_dif_EuZ)(sza)
return eps
......@@ -208,6 +301,33 @@ class calc:
def __init__(self):
pass
def PAR(self, wl, Ed):
'''
Compute instantaneous PAR from Ed spectrum.
PAR in mW m-2
PAR_quanta in µmol photon m-2 s-1
:param wl:
:param Ed:
:return:
'''
# ---------------------------------------------
# PARAMETERS
# Planck constant in J s or W s2
h = 6.6260695729e-3 # d-34
# light speed in m s-1
c = 2.99792458e0 # d8
# Avogadro Number in mol-1
Avogadro = 6.0221412927e0 # d23
hc = Avogadro * h * c
# ---------------------------------------------
idx_par = (wl >= 400) & (wl <= 700)
wl = wl[idx_par]
Ed = Ed[idx_par]
par = integrate.trapz(Ed, wl)
par_quanta = integrate.trapz(np.multiply(wl,Ed), wl) / hc
return par, par_quanta
def earth_sun_correction(self, dayofyear):
'''
Earth-Sun distance correction factor for adjustment of mean solar irradiance
......@@ -226,6 +346,40 @@ class calc:
return bidir
def angle_w(self, angle_air, n=1.334):
'''
convert above surface angle (angle_air) into sub-surface angle
:param angle_air in deg
:param n: refractive index of water
:return: sub-surface angle in deg
'''
return np.degrees(np.arcsin(np.sin(np.radians(angle_air)) / n))
def spline_2d(self, gin, arr, gout):
'''
Interpolation of a 6D array (arr) with bicubic splines on a 2D grid
corresponding to the 5th and 6th dimensions of arr.
Return 4D array interpolated on gout.
:param gin: regular 2D grid of the tabulated data (tuple/array/list of arrays)
:param arr: tabulated data (N dimensions, interpolation on N-1 and N)
:param gout: new 2D grid on which data are interpolated (with dims 2 and 3 of the same length);
(tuple/array/list of arrays)
:return: Interpolated data (1D or 3D array depending on the dimension shapes of gout
'''
N = arr.shape
interp = np.zeros(N[:-2])
for i in range(N[0]):
for j in range(N[1]):
for k in range(N[2]):
for l in range(N[3]):
interp[i, j, k, l] = interpolate.RectBivariateSpline(gin[0], gin[1], arr[i, j, k, l, ...])(
gout[0], gout[1], grid=False)
return interp
def spline_4d(self, gin, lut, gout):
'''
Interpolation with two successive bicubic splines on a regular 4D grid.
......@@ -239,7 +393,6 @@ class calc:
(tuple/array/list of arrays)
:return: Interpolated data (1D or 3D array depending on the dimension shapes of gout
'''
import scipy.interpolate as si
N = gin[0].__len__(), gin[1].__len__(), gin[2].__len__(), gin[3].__len__()
Nout = gout[0].__len__(), gout[1].__len__(), gout[2].__len__()
......@@ -247,15 +400,15 @@ class calc:
for i in range(N[0]):
for j in range(N[1]):
tmp[i, j, :] = si.RectBivariateSpline(gin[2], gin[3], lut[i, j, :, :])(gout[2], gout[3], grid=False)
tmp[i, j, :] = interpolate.RectBivariateSpline(gin[2], gin[3], lut[i, j, :, :])(gout[2], gout[3], grid=False)
if Nout[0] == Nout[1] == 1:
interp = np.ndarray(Nout[2])
for iband in range(Nout[2]):
interp[iband] = si.RectBivariateSpline(gin[0], gin[1], tmp[:, :, iband])(gout[0], gout[1], grid=False)
interp[iband] = interpolate.RectBivariateSpline(gin[0], gin[1], tmp[:, :, iband])(gout[0], gout[1], grid=False)
else:
interp = np.ndarray([Nout[0], Nout[1], Nout[2]])
for iband in range(Nout[2]):
interp[:, :, iband] = si.RectBivariateSpline(gin[0], gin[1], tmp[:, :, iband])(gout[0], gout[1],
interp[:, :, iband] = interpolate.RectBivariateSpline(gin[0], gin[1], tmp[:, :, iband])(gout[0], gout[1],
grid=True)
return interp
......@@ -13,6 +13,7 @@ from scipy.interpolate import interp1d
from utils.sunposition import sunpos
import utils.utils as u
import utils.auxdata as ua
from process.process import *
coordf = glob.glob("/DATA/OBS2CO/data/info/mesures_in_situ.csv")[0]
......@@ -23,6 +24,8 @@ awrfiles = glob.glob("/DATA/OBS2CO/data/trios/raw/aw*idpr*.csv")
# awrfiles = glob.glob("/DATA/OBS2CO/data/trios/test_setup/raw/aw*idpr*.csv")
swrfiles = glob.glob("/DATA/OBS2CO/data/trios/raw/Lu0*idpr*.csv")
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,
......@@ -36,7 +39,7 @@ idpr = '167'
# get idpr numbers
idprs = np.unique([re.findall(r'idpr(\d+)', x)[0] for x in swrfiles])
idprs = np.array(['156'])
#idprs = np.array(['170'])
# loop over idpr
for idpr in idprs:
c = coords[coords.ID_prel == int(idpr)] # .values[0]
......@@ -52,23 +55,29 @@ for idpr in idprs:
# SWR processing
# -----------------------------------------------
swr = u.swr_data(idpr, swrfiles)
if swr.file:
df, wl_swr = swr.reader(lat, lon, alt)
Rrs_swr = swr_process(df, wl_swr).process()
add_curve(ax, wl_swr, Rrs_swr.transpose().mean(axis=1), Rrs_swr.transpose().std(axis=1), label='swr', c='black')
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()
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')
# -----------------------------------------------
# AWR processing
# -----------------------------------------------
azi = 135
vza = 40
awr = u.awr_data(idpr, awrfiles)
if awr.Edf:
index_idx = [0]
azi = 135
vza = 40
d = u.data(index_idx)
Ed, wl_Ed = d.load_csv(awr.Edf)
......@@ -115,9 +124,9 @@ for idpr in idprs:
sza = sunpos(index, lat, lon, alt)[1]
df.at[index, 'sza'] = sza
rho_h = awr.get_rho_values([df.sza.min()], [vza], [azi], wl=wl)
rho15 = awr.get_rho_mobley(awr.rhoM2015, [df.sza.min()], [vza], [azi], [ws])
rho99 = awr.get_rho_mobley(awr.rhoM1999, [df.sza.min()], [vza], [azi], [ws])
rho_h = awr.get_rho_values([df.sza.mean()], [vza], [azi], wl=wl)
rho15 = awr.get_rho_mobley(awr.rhoM2015, [df.sza.mean()], [vza], [azi], [ws])
rho99 = awr.get_rho_mobley(awr.rhoM1999, [df.sza.mean()], [vza], [azi], [ws])
Rrs_h = (df.loc[:, 'Lt'] - rho_h * df.loc[:, 'Lsky']) / df.loc[:, 'Ed']
Rrs15 = (df.loc[:, 'Lt'] - rho15 * df.loc[:, 'Lsky']) / df.loc[:, 'Ed']
......@@ -125,12 +134,14 @@ for idpr in idprs:
Rrs99 = (df.loc[:, 'Lt'] - rho99 * df.loc[:, 'Lsky']) / df.loc[:, 'Ed']
# plt.figure()
add_curve(ax, wl, Rrs15.transpose().mean(axis=1), Rrs15.transpose().std(axis=1), label='M2015')
add_curve(ax, wl, Rrs99.transpose().mean(axis=1), Rrs99.transpose().std(axis=1), c='orange', label='M1999')
add_curve(ax, wl, Rrs_h.transpose().mean(axis=1), Rrs_h.transpose().std(axis=1), c='grey', label='h')
add_curve(ax, wl, Rrs15.transpose().mean(axis=1), Rrs15.transpose().std(axis=1),
label='M2015 (' + str(rho15) + ')')
add_curve(ax, wl, Rrs99.transpose().mean(axis=1), Rrs99.transpose().std(axis=1), c='orange',
label='M1999(' + str(rho99) + ')')
add_curve(ax, wl, Rrs_h.transpose().mean(axis=1), Rrs_h.transpose().std(axis=1), c='grey',
label='h(' + str(rho_h.mean()) + ')')
ax.set_title('azi=' + str(azi) + ', vza=' + str(vza)+', sza='+str(sza))
ax.set_title('azi=' + str(azi) + ', vza=' + str(vza) + ', sza=' + str(sza))
ax.legend(loc='best', frameon=False)
......@@ -138,4 +149,4 @@ for idpr in idprs:
ax.set_xlabel(r'Wavelength (nm)')
fig.savefig(os.path.join(dirfig, 'trios_awr_' + name + '_idpr' + idpr + '.png'))
plt.close()
Lt.index.names
......@@ -71,7 +71,7 @@ fig.subplots_adjust(left=0.1, right=0.9, hspace=.5, wspace=0.65)
i=0
for azi, Lt1 in Lt0.groupby(level=2):
for vza,Lt in Lt1.groupby(level=1):
ax=axs.flat[i]
ax = axs.flat[i]
i=i+1
print(azi,vza)
......
......@@ -10,7 +10,7 @@ class awr_data:
Above-water radiometry
'''
def __init__(self, idpr, files):
def __init__(self, idpr=None, files=None):
# ''' get file names for Ed, Lsky and Lt data'''
self.file = list(filter(lambda x: 'idpr' + idpr in x, files))
file = self.file
......@@ -41,7 +41,7 @@ class awr_data:
df = pd.DataFrame()
# ''' read files with pandas format '''
d=data(index_idx)
d = data(index_idx)
Ed, wl_Ed = d.load_csv(self.Edf)
Lsky, wl_Lsky = d.load_csv(self.Lskyf)
Lt, wl_Lt = d.load_csv(self.Ltf)
......@@ -90,7 +90,7 @@ class iwr_data:
self.Luzf = list(filter(lambda x: 'Luz' in x, file))
self.idpr = idpr
def reader(self, lat=None, lon=None, alt=0, name=''):
def reader(self, lat, lon, alt=0, name='', delta_Lu_depth=0, delta_Edz_depth=0):
'''
Read above-water data files for a given acquisition series (idpr),
merge the different data types:
......@@ -105,6 +105,9 @@ class iwr_data:
:param lat: latitude (decimal)
:param lon: longitude (decimal)
:param alt: altitude (m)
:param delta_Lu_depth: adjustment of actual depth for Lu sensor (distance from depth sensor);
in meters for depth counted positively
:param delta_Edz_depth: similar to delta_Lu_depth for Edz sensor
:param idpr: ID of the acquisition series
:return:
'''
......@@ -112,29 +115,42 @@ class iwr_data:
df = pd.DataFrame()
# ''' read files with pandas format '''
Ed, wl_Ed = self.load_csv(self.Edf)
Edz, wl_Edz = self.load_csv(self.Edzf)
Luz, wl_Luz = self.load_csv(self.Luzf)
d = data([1, 0])
Ed, wl_Ed = d.load_csv(self.Edf)
Edz, wl_Edz = d.load_csv(self.Edzf)
Luz, wl_Luz = d.load_csv(self.Luzf)
#mask negative values TODO save number of discarded data
# Ed.mask(Ed<0,inplace=True)
# Edz.mask(Edz<0,inplace=True)
# Luz.mask(Luz<0,inplace=True)
# copy depth data to Ed frame on date index
Ed.index = Ed.index.droplevel(level=1)
# Ed.index = Ed.index.droplevel(level=1)
# ''' interpolate Ed and Lsky data upon Lt wavelength'''
#''' interpolate Ed and Lsky data upon Lt wavelength'''
wl = wl_Luz
Luz.columns = pd.MultiIndex.from_tuples(zip(['Luz'] * len(wl), wl), names=['param', 'wl'])
Luz.columns = pd.MultiIndex.from_tuples(list(zip(['Luz'] * len(wl), wl)), names=['param', 'wl'])
intEd = interp1d(wl_Ed, Ed.values, fill_value='extrapolate')(wl)
newEd = pd.DataFrame(index=Ed.index,
columns=pd.MultiIndex.from_tuples(zip(['Ed'] * len(wl), wl), names=['param', 'wl']),
newEd = pd.DataFrame(index=Ed.index.get_level_values(0),
columns=pd.MultiIndex.from_tuples(list(zip(['Ed'] * len(wl), wl)), names=['param', 'wl']),
data=intEd)
intEdz = interp1d(wl_Edz, Edz.values, fill_value='extrapolate')(wl)
newEdz = pd.DataFrame(index=Edz.index, columns=pd.MultiIndex.from_tuples(zip(['Edz'] * len(wl), wl),
newEdz = pd.DataFrame(index=Edz.index, columns=pd.MultiIndex.from_tuples(list(zip(['Edz'] * len(wl), wl)),
names=['param', 'wl']), data=intEdz)
# merge sensor data on time
# correct depth data for sensor to sensor distance
Luz.reset_index(level=1, inplace=True)
Luz.iloc[:, 0] = Luz.iloc[:, 0] + delta_Lu_depth
# newEd.reset_index(level=1,inplace=True)
newEdz.reset_index(level=1, inplace=True)
newEdz.iloc[:, 0] = newEdz.iloc[:, 0] + delta_Edz_depth
# merge sensor data on time
df = pd.merge_asof(Luz, newEd, left_index=True, right_index=True, tolerance=pd.Timedelta("2 seconds"),
direction="nearest")
df = pd.merge_asof(df, newEdz, left_index=True, right_index=True, suffixes=('_Luz', '_Edz'),
......@@ -150,26 +166,26 @@ class iwr_data:
df.at[index, 'sza'] = sza
df['idpr', ''] = self.idpr
df['name', ''] = self.name
df['name', ''] = name
return df, wl
def load_csv(self, file):
dateparse = lambda x: pd.datetime.strptime(x, '%Y-%m-%d %H:%M:%S')
if len(file) > 1:
print('Warning! Multiple files found but only one expected, process first file of the list:')
print(file)
file = file[0]
df = pd.read_csv(file, sep=';', index_col=[1, 0], na_values=['-NAN'])
df = df.dropna(axis=1, how='all').dropna(axis=0, how='all')
df.index = df.index.set_levels([pd.to_datetime(df.index.levels[0]), df.index.levels[1]])
df.columns = df.columns.astype('float') # str.extract('(\d+)',expand=False).astype('float')
# resort to get data in increasing time order
df.sort_index(inplace=True, level=0)
wl = df.columns
return df, wl
# def load_csv(self, file):
#
# dateparse = lambda x: pd.datetime.strptime(x, '%Y-%m-%d %H:%M:%S')
# if len(file) > 1:
# print('Warning! Multiple files found but only one expected, process first file of the list:')
# print(file)
# file = file[0]
# df = pd.read_csv(file, sep=';', index_col=[1, 0], na_values=['-NAN'])
# df = df.dropna(axis=1, how='all').dropna(axis=0, how='all')
# df.index = df.index.set_levels([pd.to_datetime(df.index.levels[0]), df.index.levels[1]])
# df.columns = df.columns.astype('float') # str.extract('(\d+)',expand=False).astype('float')
# # resort to get data in increasing time order
# df.sort_index(inplace=True, level=0)
# wl = df.columns
#
# return df, wl
class swr_data:
......@@ -235,28 +251,28 @@ class swr_data:
class data:
def __init__(self,index_idx=[0]):
def __init__(self, index_idx=[0]):
# first position should be datetime index
# followed by the other parameters used for indexing (e.g. azimuth, view angle)
self.index_idx=index_idx
self.index_idx = index_idx
pass
def load_csv(self, file):