commit ceac9f1756b2ba970dc434d7c25412db94bf041d
Author: tangchao <735056338@qq.com>
Date:   Wed May 25 23:32:32 2022 +0800

    第一次提交

diff --git a/.gitignore b/.gitignore
new file mode 100644
index 0000000..09b3c4a
--- /dev/null
+++ b/.gitignore
@@ -0,0 +1,164 @@
+# 唐超添加
+/.idea
+
+
+# Byte-compiled / optimized / DLL files
+__pycache__/
+*.py[cod]
+*$py.class
+
+# C extensions
+*.so
+
+# Distribution / packaging
+.Python
+build/
+develop-eggs/
+dist/
+downloads/
+eggs/
+.eggs/
+lib/
+lib64/
+parts/
+sdist/
+var/
+wheels/
+share/python-wheels/
+*.egg-info/
+.installed.cfg
+*.egg
+MANIFEST
+
+# PyInstaller
+#  Usually these files are written by a python script from a template
+#  before PyInstaller builds the exe, so as to inject date/other infos into it.
+*.manifest
+*.spec
+
+# Installer logs
+pip-log.txt
+pip-delete-this-directory.txt
+
+# Unit test / coverage reports
+htmlcov/
+.tox/
+.nox/
+.coverage
+.coverage.*
+.cache
+nosetests.xml
+coverage.xml
+*.cover
+*.py,cover
+.hypothesis/
+.pytest_cache/
+cover/
+
+# Translations
+*.mo
+*.pot
+
+# Django stuff:
+*.log
+local_settings.py
+db.sqlite3
+db.sqlite3-journal
+
+# Flask stuff:
+instance/
+.webassets-cache
+
+# Scrapy stuff:
+.scrapy
+
+# Sphinx documentation
+docs/_build/
+
+# PyBuilder
+.pybuilder/
+target/
+
+# Jupyter Notebook
+.ipynb_checkpoints
+
+# IPython
+profile_default/
+ipython_config.py
+
+# pyenv
+#   For a library or package, you might want to ignore these files since the code is
+#   intended to run in multiple environments; otherwise, check them in:
+# .python-version
+
+# pipenv
+#   According to pypa/pipenv#598, it is recommended to include Pipfile.lock in version control.
+#   However, in case of collaboration, if having platform-specific dependencies or dependencies
+#   having no cross-platform support, pipenv may install dependencies that don't work, or not
+#   install all needed dependencies.
+#Pipfile.lock
+
+# poetry
+#   Similar to Pipfile.lock, it is generally recommended to include poetry.lock in version control.
+#   This is especially recommended for binary packages to ensure reproducibility, and is more
+#   commonly ignored for libraries.
+#   https://python-poetry.org/docs/basic-usage/#commit-your-poetrylock-file-to-version-control
+#poetry.lock
+
+# pdm
+#   Similar to Pipfile.lock, it is generally recommended to include pdm.lock in version control.
+#pdm.lock
+#   pdm stores project-wide configurations in .pdm.toml, but it is recommended to not include it
+#   in version control.
+#   https://pdm.fming.dev/#use-with-ide
+.pdm.toml
+
+# PEP 582; used by e.g. github.com/David-OConnor/pyflow and github.com/pdm-project/pdm
+__pypackages__/
+
+# Celery stuff
+celerybeat-schedule
+celerybeat.pid
+
+# SageMath parsed files
+*.sage.py
+
+# Environments
+.env
+.venv
+env/
+venv/
+ENV/
+env.bak/
+venv.bak/
+
+# Spyder project settings
+.spyderproject
+.spyproject
+
+# Rope project settings
+.ropeproject
+
+# mkdocs documentation
+/site
+
+# mypy
+.mypy_cache/
+.dmypy.json
+dmypy.json
+
+# Pyre type checker
+.pyre/
+
+# pytype static type analyzer
+.pytype/
+
+# Cython debug symbols
+cython_debug/
+
+# PyCharm
+#  JetBrains specific template is maintained in a separate JetBrains.gitignore that can
+#  be found at https://github.com/github/gitignore/blob/main/Global/JetBrains.gitignore
+#  and can be added to the global gitignore or merged into this file.  For a more nuclear
+#  option (not recommended) you can uncomment the following to ignore the entire idea folder.
+#.idea/
\ No newline at end of file
diff --git a/Setup_Window.m b/Setup_Window.m
new file mode 100644
index 0000000..8d51066
--- /dev/null
+++ b/Setup_Window.m
@@ -0,0 +1,72 @@
+function [windows] = Setup_Window()
+    
+    windows = [];
+    % -----Attention: fitting window range for sFLD��3FLD��SFM��SVD��DOAS;----- %
+    % ����sFLD, iFLD, pFLD
+    FLD = [];
+    FLD.outerH = [656.0, 656.2];     % �˵㲨�Σ�ȡ��ֵ
+    FLD.outerB = [687.0, 687.2];   
+    FLD.outerW = [717.0, 717.2];  
+    FLD.outerA = [758.6, 758.8]; 
+    FLD.windowH = [656, 659];         % ��Χ������Сֵ
+    FLD.windowB = [686, 689];    
+    FLD.windowW = [717, 722]; 
+    FLD.windowA = [758, 770]; 
+    windows.FLD = FLD;
+
+    % 3FLD
+    FLD3 = [];
+    FLD3.FLD3_leftH  = [656.0, 656.2];    FLD3.FLD3_rightH = [659.0, 659.2];  
+    FLD3.FLD3_leftB  = [687.0, 687.2];    FLD3.FLD3_rightB = [689.0, 689.2];    
+    FLD3.FLD3_leftW  = [717.0, 717.2];    FLD3.FLD3_rightW = [722.0, 722.2]; 
+    FLD3.FLD3_leftA  = [758.6, 758.8];    FLD3.FLD3_rightA = [770.0, 770.2];
+    FLD3.windowH = [656, 659];           % ��Χ������Сֵ
+    FLD3.windowB = [686, 689];    
+    FLD3.windowW = [717, 722]; 
+    FLD3.windowA = [758, 770]; 
+    windows.FLD3 = FLD3;
+    
+    % SFM
+    SFM = [];
+    SFM.SFM_windowH = [656, 667];       % ���η�Χ
+    SFM.SFM_windowB = [687, 698];
+    SFM.SFM_windowW = [717, 728];
+    SFM.SFM_windowA = [759, 770];
+    SFM.SFM_Ha =  [656.0 656.2];       % ij����ӫ��ֵ,SFM��SVD�㷨ֱ�Ӷ�����
+    SFM.SFM_O2B = [687.0 687.2];
+    SFM.SFM_H2O = [720.0 720.2];
+    SFM.SFM_O2A = [760.0 760.2];
+    windows.SFM = SFM;
+    % SVD
+    SVD = [];
+    SVD.SVD_windowH  = [655, 670];     
+    SVD.SVD_windowB  = [675, 710];      % wide O2B and O2A fitting window; % e.g. narrow O2A fitting window:[759.86 762.79];
+    SVD.SVD_windowW  = [700, 740];
+    SVD.SVD_windowA  = [740, 780];  
+    SVD.SVD_Ha =  [656.0 656.2];       % ij����ӫ��ֵ,SFM��SVD�㷨ֱ�Ӷ�����
+    SVD.SVD_O2B = [687.0 687.2];
+    SVD.SVD_H2O = [720.0 720.2];
+    SVD.SVD_O2A = [760.0 760.2];
+    % Example SVD settings
+    % if rolling > 0, abs(rolling) = the number of spectra in a centered moving window used for training the SVD
+    % if rolling = 0, it will use all spectra from the day for training the SVD
+    % if rolling < 0, it will use that number of spectra but evenly dispersed throughout the day for training the SVD
+    SVD.rolling = 0;                    % moving window of 5 spectra 
+    windows.SVD = SVD;
+    
+    % DOAS, F_SFM, ���ɷַ�����Ҳ����ȡ���η�Χ�ڷ�������Сֵ       
+    Others = [];                           % ���ַ��������õIJ��η�Χ��ȡ��Сֵ����� 
+    Others.WindowH = [656, 659];           % Attention: ע���޸�
+    Others.WindowB = [686, 689];    
+    Others.WindowW = [717, 722]; 
+    Others.WindowA = [758, 770]; 
+    windows.Others  = Others;
+    
+    % -----Attention: Modify the range of spectrum;----- %
+    % --�õ�������������նȣ������Ⱥͷ�����
+    windows.wl_680 = [676, 680];    % 676-680nm, red
+    windows.wl_705 = [703, 707];    % 703-707nm, red
+    windows.wl_740 = [736, 740];    % 736-740nm, rededge
+    windows.wl_770 = [770, 774];    % 770-774nm, near_red
+    windows.wl_800 = [795, 799];    % 795-799nm, near_red
+end
\ No newline at end of file
diff --git a/sif_methods.py b/sif_methods.py
new file mode 100644
index 0000000..215ade8
--- /dev/null
+++ b/sif_methods.py
@@ -0,0 +1,310 @@
+import numpy as np
+from scipy.special import eval_legendre
+from scipy import optimize
+# 所有提取方法均基于xarray数据，xarray数据的变量与波长、时间、光谱仪绑定
+def cal_inside_bands_ave(data):
+    '''
+    根据多个光谱找出窗口内数据最低点对应波长
+    '''
+    sky_spec = data.sky
+    veg_spec = data.veg
+    wvl_inside_band_l = np.mean(veg_spec.idxmin(dim='Wavelength')).values
+    wvl_inside_band_e = np.mean(sky_spec.idxmin(dim='Wavelength')).values
+    return[wvl_inside_band_l,wvl_inside_band_e]
+
+def cal_outside_values_mean(data,outer):
+    '''
+    计算肩部窗口的均值
+    '''
+    _data = data.where((data.Wavelength>outer[0])&(data.Wavelength<outer[1]),drop=True)
+    wvl_outer_mean = _data['Wavelength'].mean().values
+    _mean = _data.mean(dim = 'Wavelength')
+    E_out = _mean.sky.values
+    L_out = _mean.veg.values
+    return L_out,E_out,wvl_outer_mean
+
+def sfld(data,wl_range,outer):
+
+    """
+    Standard FLD (Huaize Feng)
+
+    input:
+          data,xarray dataset
+          wl_range, a window around Fraunhofer lines position, [start,end]
+          outer, a window outside the Fraunhofer lines, [start,end]
+    """
+    sif = []
+    refl = []
+    nmeas_ = data.Measures.size
+    
+    data = data.where((data.Wavelength>wl_range[0])&(data.Wavelength<wl_range[1]),drop=True)
+    [wvl_inside_band_l,wvl_inside_band_e]=cal_inside_bands_ave(data)
+    for i in range(0,nmeas_):
+        _data = data.isel(Measures=i)
+        veg_out,sky_out,_ = cal_outside_values_mean(_data,outer)
+        veg_in = _data.veg.sel(Wavelength = wvl_inside_band_l,method='nearest').values
+        sky_in = _data.sky.sel(Wavelength = wvl_inside_band_e,method='nearest').values
+        _sif = (sky_out*veg_in - sky_in*veg_out)/(sky_out - sky_in) 
+        _refl = (veg_in - _sif)/ sky_in
+        sif.append(_sif)
+        refl.append(_refl)
+    return[sif,refl]
+
+def fld3(data,wl_range,outer_left,outer_right):
+    """
+    3FLD (Huaize Feng)
+        input:
+          data,xarray dataset
+          wl_range, a window around Fraunhofer lines position, [start,end]
+          outer_*, a window outside the Fraunhofer lines, [start,end]
+              outer_left, left window outside the Fraunhofer lines
+              outer_right, left window outside the Fraunhofer lines
+        return:
+          list[sif,reflectance]
+    """
+    sif = []
+    refl = []
+    nmeas_ = data.Measures.size
+    
+    data = data.where((data.Wavelength>wl_range[0])&(data.Wavelength<wl_range[1]),drop=True)
+    [wvl_inside_band_l,wvl_inside_band_e]=cal_inside_bands_ave(data)
+    for i in range(0,nmeas_):
+        _data = data.isel(Measures=i)
+        veg_out_left,sky_out_left,wvl_outer_left = cal_outside_values_mean(_data,outer_left)
+        veg_out_right,sky_out_right,wvl_outer_right = cal_outside_values_mean(_data,outer_right)
+        veg_in = _data.veg.sel(Wavelength = wvl_inside_band_l,method='nearest').values
+        sky_in = _data.sky.sel(Wavelength = wvl_inside_band_e,method='nearest').values
+        
+        # 根据离吸收峰的距离的反比进行赋权
+        wight_left = (wvl_outer_right - wvl_inside_band_e)/(wvl_outer_right - wvl_outer_left)
+        wight_right = (wvl_inside_band_e - wvl_outer_left)/(wvl_outer_right - wvl_outer_left)
+        
+        _sif =  (veg_in - (sky_in/((wight_left*sky_out_left) + (wight_right*sky_out_right))) * ((wight_left*veg_out_left) + (wight_right*veg_out_right))) / (1-(sky_in / ((wight_left*sky_out_left) + (wight_right*sky_out_right))))
+        _refl = (veg_in - _sif)/ sky_in
+        sif.append(_sif)
+        refl.append(_refl)
+    return[sif,refl]        
+
+def sfm(data,wl_range,band):
+
+    """
+    Spectral Fitting Method (Huaize Feng)
+    input:
+          data,
+              xarray dataset
+          wl_range,[start,end]
+              a window around Fraunhofer lines position, ,about 10 nm
+              for fitting a ployniam or gussian function.
+          band,float/int
+              exact position of the Fraunhofer line
+    return:
+          list[sif,reflectance,rmse,B]
+              sif, numpy array, sif at the Fraunhofer line for each measurement
+              reflectance, numpy array
+              rmse, RMSE
+              B, numpy array, the parameters of the fitting equation
+                  [a,b,c,d,e,f]: a, b, c for polynimal refelectance
+                                 d, float, MAX sif value [0,10] mw/m2/nm/sr
+                                 e, position of MAX sif value - wl_range, [0,wavelength.size], nm
+                                 f, full width at half maximum, [0,wavelength.size*5]
+    """
+    sif,refl,rmse,B = [],[],[],[]
+    data = data.where((data.Wavelength>wl_range[0])&(data.Wavelength<wl_range[1]),drop=True)
+    abosorb_line_position = np.where(data.Wavelength == data.Wavelength.sel(Wavelength=760,method='nearest').values)
+    _nmeas = data.Measures.size
+    _nwvl = data.Wavelength.size
+    _x = (data.Wavelength -np.min(data.Wavelength)).values
+    poly_refl = [_x**2,_x,np.ones(_nwvl)]
+    poly_sif = [_x**2,_x,np.ones(_nwvl)]
+   
+    for i in range(0,_nmeas):
+        sky_spec = data.sky[i].values
+        veg_spec = data.veg[i].values
+        _X = np.concatenate([poly_refl*sky_spec,poly_refl]).T
+        _B = np.linalg.lstsq(_X,veg_spec,rcond=-1)
+        _sif = np.array(poly_sif).T.dot(_B[0][-3:])
+        _refl = (veg_spec - _sif) / (sky_spec+0.000001)
+        sif.append(_sif[abosorb_line_position][0])
+        refl.append(_refl[abosorb_line_position][0])
+        rmse.append(np.sqrt(np.sum(_B[1]**2)/_nwvl))
+        B.append(_B[0])
+    return [sif,refl,rmse,B]
+
+
+def f(x, a, b, c, d, e, f):
+    sky_spec = x['sky_spec']
+    _x = x['_x']
+    refl = a * _x**2*sky_spec + b* _x*sky_spec + c*sky_spec
+    sif = d*np.exp(-(_x-e)**2/f)
+    y_hat = refl + sif
+    return y_hat
+
+def f_cal(X, a,b,c,d,e,f):
+    sky_spec = X['sky_spec']
+    _x = X['_x']
+    refl = a * _x**2*sky_spec + b* _x*sky_spec + c*sky_spec
+    sif = d*np.exp(-(_x-e)**2/f)
+    y_hat = refl + sif
+    return sif,refl
+
+def sfm_gaussian(data,wl_range,band=760):
+    
+    '''
+    Spectral Fitting Method (Gaussian Ver. (Huaize Feng))
+    input:
+          data,
+              xarray dataset
+          wl_range,[start,end]
+              a window around Fraunhofer lines position, ,about 10 nm
+              for fitting a polynomial or gaussian function.
+          band,float/int
+              exact position of the Fraunhofer line
+    return:
+          list[sif,reflectance,rmse,B]
+              sif, numpy array, sif at the Fraunhofer line for each measurement
+              reflectance, numpy array
+              rmse, RMSE
+              B, numpy array, the parameters of the fitting equation
+                  [a,b,c,d,e,f]: a, b, c for polynimal refelectance
+                                 d, float, MAX sif value [0,10] mw/m2/nm/sr
+                                 e, position of MAX sif value - wl_range, [0,wavelength.size], nm
+                                 f, full width at half maximum
+    '''
+    sif,refl,rmse,B = [],[],[],[]
+    data = data.where((data.Wavelength>wl_range[0])&(data.Wavelength<wl_range[1]),drop=True)
+    abosorb_line_position = np.where(data.Wavelength == data.Wavelength.sel(Wavelength=band,method='nearest').values)
+    _nmeas = data.Measures.size
+    _nwvl = data.Wavelength.size
+    _x = (data.Wavelength -np.min(data.Wavelength)).values
+    poly_refl = [_x**2,_x,np.ones(_nwvl)]
+    poly_sif = [_x**2,_x,np.ones(_nwvl)]
+    
+    bounds=((-np.inf, -np.inf,-np.inf, 0, 0, 0), (np.inf, np.inf,np.inf, 10, _nwvl, _nwvl*5))
+    
+    for i in range(0,_nmeas):
+        sky_spec = data.sky[i].values
+        veg_spec = data.veg[i].values
+        _X = {'sky_spec':sky_spec,'_x':_x}
+        _B = optimize.curve_fit(f, _X, veg_spec,bounds=bounds)[0]
+#         print(_B)
+        
+        _sif,_ = f_cal(_X,*_B)
+#         print(_sif)
+        _rmse = np.sqrt(np.sum((f(_X,*_B)-veg_spec)**2)/_nwvl)
+        _refl = (veg_spec - _sif) / sky_spec
+        sif.append(_sif[abosorb_line_position][0])
+        refl.append(_refl[abosorb_line_position][0])
+        rmse.append( _rmse)
+        B.append(_B[0])
+    return [sif,refl,rmse,B]
+
+def doas(data,wl_range,band=760):
+
+    """
+    Spectral Fitting Method (Huaize Feng)
+    input:
+          data,
+              xarray dataset
+          wl_range,[start,end]
+              a window around Fraunhofer lines position, ,about 10 nm
+              for fitting a ployniam or gussian function.
+          band,float/int
+              exact position of the Fraunhofer line
+    return:
+          list[sif,reflectance,rmse,B]
+              sif, numpy array, sif at the Fraunhofer line for each measurement
+              reflectance, numpy array
+              rmse, RMSE
+              B, numpy array, the parameters of the fitting equation
+                  [a,b,c,d,e,f]: a, b, c for polynimal refelectance
+                                 d, float, MAX sif value [0,10] mw/m2/nm/sr
+                                 e, position of MAX sif value - wl_range, [0,wavelength.size], nm
+                                 f, full width at half maximum, [0,wavelength.size*5]
+    """
+    sif,refl,rmse,B = [],[],[],[]
+    data = data.where((data.Wavelength>wl_range[0])&(data.Wavelength<wl_range[1]),drop=True)
+    absorb_line_position = np.where(data.Wavelength == data.Wavelength.sel(Wavelength=band,method='nearest').values)
+    _nmeas = data.Measures.size
+    _wvl = data.Wavelength.values
+    _nwvl = _wvl.size
+    _x = np.interp(_wvl, (_wvl.min(),_wvl.max()),(-1,1))
+    # normalize standard SIF template by the mean value
+    _hf = data.hf.values/(np.mean(data.hf.values))
+    # create base function by legendre polynomial equations
+    _ref_base = np.array([eval_legendre(n, _x) for n in np.arange(6)])
+    
+    for i in range(0,_nmeas):
+        print(data.Measures[i].values)
+        sky_spec = data.sky[i].values
+        veg_spec = data.veg[i].values
+        _X = np.concatenate([_ref_base,(_hf/(veg_spec+0.000001111)).reshape(1,-1)]).T
+        _y = np.log(veg_spec+0.000001111) - np.log(sky_spec+0.00000111111)
+        _B = np.linalg.lstsq(_X,_y,rcond=-1)
+        _sif = _hf * _B[0][-1]
+        _refl = (veg_spec - _sif) / (sky_spec+0.0000011111)
+        sif.append(_sif[absorb_line_position][0])
+        refl.append(_refl[absorb_line_position][0])
+        rmse.append(np.sqrt(np.sum(_B[1]**2)/_nwvl))
+        B.append(_B[0])
+        
+    return [sif,refl,rmse,B]
+
+def svd(data,wl_range,band=760,num_vector=20,pow_of_refl = 5):
+
+    """
+    Singular Vector Decomposition (Huaize Feng)
+    input:
+          data,
+              xarray dataset
+          wl_range,[start,end]
+              a window around Fraunhofer lines position, ,about 10 nm
+              for fitting a polynomial or gaussian function.
+          band,float/int
+              exact position of the Fraunhofer line
+          num_vector,int
+                pca保留的观测维度（个数），相当于去除大气的噪音
+          pow_of_refl,int
+                反射率所用多项式的最高次数
+              Highest power of the polynomial equation denoting reflectance
+    return:
+          list[sif,reflectance,rmse,B]
+              sif, numpy array, sif at the Fraunhofer line for each measurement
+              reflectance, numpy array
+              rmse, RMSE
+              B, numpy array, the parameters of the fitting equation
+                  [a,b,c,d,e,f]: a, b, c for polynimal refelectance
+                                 d, float, MAX sif value [0,10] mw/m2/nm/sr
+                                 e, position of MAX sif value - wl_range, [0,wavelength.size], nm
+                                 f, full width at half maximum, [0,wavelength.size*5]
+    """
+    sif,refl,rmse,B = [],[],[],[]
+    data = data.where((data.Wavelength>wl_range[0])&(data.Wavelength<wl_range[1]),drop=True)
+    abosorb_line_position = np.where(data.Wavelength == data.Wavelength.sel(Wavelength=band,method='nearest').values)
+    _nmeas = data.Measures.size
+    _wvl = data.Wavelength.values
+    _nwvl = _wvl.size
+    _hf = data.hf.values  # 重采样后的标准sif
+    
+    u, s, vh = np.linalg.svd(data.sky)  # Svd分解----------------------------------------------------------------------------------------------------
+    v = -vh
+
+    tmp1=(_wvl-np.mean(_wvl)).reshape(-1,1)
+    tmp2=np.arange(pow_of_refl+1)
+    p1 = np.power(tmp1, tmp2)*v[0].reshape(-1,1)
+    p2 = np.power((_wvl-np.mean(_wvl)).reshape(-1,1), np.arange(pow_of_refl+1))*v[1].reshape(-1,1)
+    p3 = v[:,2:num_vector]
+    _X = np.concatenate([p1,p2,p3,_hf.reshape(-1,1)],axis=1)
+    
+    for i in range(0,_nmeas):
+        sky_spec = data.sky[i].values
+        veg_spec = data.veg[i].values
+        _y = veg_spec
+        _B = np.linalg.lstsq(_X,_y,rcond=-1)  # https://www.zhihu.com/question/40540185?sort=created
+        _sif = _hf * _B[0][-1]
+        _refl = (veg_spec - _sif) / (sky_spec + 0.0000001)
+        sif.append(_sif[abosorb_line_position][0])
+        refl.append(_refl[abosorb_line_position][0])
+        rmse.append(np.sqrt(np.sum(_B[1]**2)/_nwvl))
+        B.append(_B[0])
+    return [sif,refl,rmse,B]
+
diff --git a/sif_retrieval.py b/sif_retrieval.py
new file mode 100644
index 0000000..a21956b
--- /dev/null
+++ b/sif_retrieval.py
@@ -0,0 +1,251 @@
+# Multi Processing Version.
+# 2022/1/14 14:52
+
+# 唐超：如何将此文件打包成一个exe
+#  （1）cd D:\csharp_vs2017\easySif\sifAlgorithm_Python_exe
+#  （2）D:\software\Anaconda3\envs\python_36\Scripts\pyi-makespec.exe -F D:\PycharmProjects\sif\sif_retrieval.py
+#  （3）在生成的spec文件中加入import sys
+#                            sys.setrecursionlimit(5000)
+#  （4）D:\software\Anaconda3\envs\python_36\Scripts\pyinstaller.exe  D:\csharp_vs2017\easySif\sifAlgorithm_Python_exe\sif_retrieval.spec
+
+import os
+import sys
+from glob import glob
+import ast
+import csv
+import sif_methods
+import argparse
+from multiprocessing import Pool
+import multiprocessing
+from scipy.interpolate import splev, splrep
+import numpy as np
+import pandas as pd
+
+try:
+    import xarray as xr
+except Exception as e:
+    print(e)
+    cmd = 'pip install xarray -i https://pypi.tuna.tsinghua.edu.cn/simple --user'
+    os.system(cmd)
+    import xarray as xr
+finally:
+    pass
+
+
+def intep(standard, x2):
+    '''
+    插值，根据给定的波长序列将标准SIF重新插值
+    输入：标准SIF，两列，一列波长一列值；给定新的波长序列
+    输出：插值后的SIF值
+    '''
+    spl = splrep(standard.index, standard)
+    y2 = splev(x2, spl)
+    return y2
+
+
+def sel_method(x):
+    '''
+    选择方法，
+    输入：方法字符串
+    输出：方法
+    '''
+    return {
+        'sfld': sif_methods.sfld,
+        '3fld': sif_methods.fld3,
+        'svd': sif_methods.svd,
+        'sfm': sif_methods.sfm,
+        'sfm_gaussian': sif_methods.sfm_gaussian,
+        'doas': sif_methods.doas
+    }.get(x, sif_methods.sfld)
+
+
+def read_files(f):
+    '''
+    使用csv.reader读取文件
+    '''
+    with open(f) as csvfile:
+        data = list(csv.reader(csvfile))  # 读csv并转化为list
+#         data = np.array(data)
+        data = np.array(data, dtype=object)  # 转化为np数组
+        data = [np.array(d) for d in data]
+        csvfile.close()
+    return data
+
+
+def Convert(lst):
+    '''
+    针对<键名 键值...>组织形式的文件头解析
+    '''
+    # tc
+    #     res_dct = {lst[i]: lst[i + 1] for i in range(0, lst.size-1, 2)}
+    res_dct = {lst[i]: lst[i + 1] for i in range(0, len(lst)-1, 2)}
+    return res_dct
+
+# parsing
+
+
+def data_parser(data, dt):
+    shp = np.arange(len(data)).tolist()
+    _iter = iter(shp)
+
+    # 读取csv的元信息（第一行）
+    data[next(_iter)][1:]  # 略过csv的第一行
+    header_dict = {'Measure': dt}
+    # header_dict = Convert(data[next(_iter)][1:])
+    header_dict.update(Convert(data[next(_iter)]))
+
+    # 读取波长等信息
+    ls = []
+    for i in range(int(header_dict['TotalSpectrometer'])):
+        info = Convert(data[next(_iter)][1:])  # 略过第3行
+        wvl = np.array(data[next(_iter)][1:])
+        info['Wavelength'] = (wvl[wvl != '']).astype('float')
+        info.update(header_dict)
+#         info = my_dict
+        ls.append(info)
+
+    # 读取光谱数据?????????????????????????????????????????????????????????????????????????????????????????????????????????????????????
+    data = data[next(_iter)+1:]
+    data = pd.DataFrame(data)
+    data[['spec', 'point']] = data[0].str.split('_', expand=True)  #??
+    data = data.set_index(['spec', 'point'])
+    data = data.iloc[:, 3:]
+    data.columns.name = 'Wavelength'
+    data = data.astype('float')
+    data = data.T.unstack().to_xarray()
+    # _ ['spec'] = f
+    data['Measures'] = dt
+    data = data.expand_dims('Measures')
+    return [data, ls]  # data是光谱数据，ls是元信息
+
+
+def map_files(f):
+    '''
+    遍历文件
+    '''
+    dt = pd.to_datetime(f.split('\\')[-1][:-4], format='%Y_%m_%d_%H_%M_%S')
+    data = read_files(f)
+    data = data_parser(data, dt)
+    #     par_ls.append(par)
+    return data    # data是列表：第一个是光谱数据，第二个是元信息
+
+# def map_files(i):
+#     result = i * i
+#     return result
+
+
+def processing(standard_sif, folder, out_file, pars, data, header, sky_p='P1', method='sfld'):
+    '''
+    提取算法调用
+    '''
+    # 参数解析为list
+    pars = ast.literal_eval(pars)
+    standard_sif = pd.read_csv(standard_sif, index_col=[0])
+    point_ls = {}
+    ls = []
+    # 遍历光谱仪
+    for spec in data['spec'].values:
+        if "SIF" not in spec:
+            continue
+        else:
+            # 
+            _ = data.sel(spec=spec)  # ???????????????
+            # 根据光谱仪提取对应的波长
+            wvl = header['Wavelength'][spec]
+            size = wvl.size
+            # 使用波长的长度对数据进行截取（减少多光谱仪数据不一致导致的空值）
+            _ = _.isel(Wavelength=xr.DataArray(
+                np.arange(0, size), dims="Wavelength"))
+            _['Wavelength'] = wvl
+        #     _ = _.where((_.Wavelength>731.3)&(_.Wavelength<782),drop=True)
+
+            sky = _.sel(point=sky_p, drop=True).rename('sky')
+            for p in _.point:
+                if p == sky_p:
+                    continue
+                else:
+                    veg = _.sel(point=p, drop=True).rename('veg')
+                    input_each = xr.merge([sky, veg])
+                    _hf = intep(standard_sif, input_each.Wavelength.values)  # 将标准sif插值匹配到数据的波长
+                    input_each['hf'] = (['Wavelength'], _hf)
+                    # 调用方法
+                    retr_method = sel_method(method)
+                    print('Running {} method on {} of spectrometer {}'.format(
+                        method, p.values, spec), flush=True)
+                    print('Processing ...', flush=True)
+                    sif = retr_method(input_each, *pars)[0]#---------------------------------------------------------------------------------------------
+                    point_ls.update({str(p.values): sif})
+            _sif = pd.DataFrame(point_ls)
+            point_ls = {}
+            _sif.index = sky.Measures.values
+            _sif.index.name = 'Measures'
+            _sif = _sif.to_xarray()
+            _sif = _sif.assign_coords(spec=spec)
+            ls.append(_sif)
+    sif = xr.merge(ls)
+    sif.to_dataframe().to_csv(out_file)
+
+
+
+
+# # 需要输入的参数
+# # 标准SIF曲线
+# standard_sif = 'standard_sif.csv'
+# # 天空辐射对应的点
+# sky_p = 'P1'
+# # 存放文件夹
+# folder = '2021_12_21'
+# # 提取方法需要的输入参数，用括号括起来
+# pars = ([740,770],760)
+
+# processing(standard_sif,folder,outfile,pars,sky_p='P1',method='svd')
+
+def parse_args():
+    parser = argparse.ArgumentParser(description="SIF retrieval process")
+    parser.add_argument("standard_sif")
+    parser.add_argument("folder")
+    parser.add_argument("outfile")
+    parser.add_argument("pars")
+    parser.add_argument("sky_p", default='P1')
+    parser.add_argument("method", default='sfld')
+    args = parser.parse_args()
+    return args
+
+
+def main():
+    # 解析命令行参数
+    inputs = parse_args()
+
+    # 读取数据
+    cpus = os.cpu_count()
+    pool = Pool(cpus)
+    # print('tc-----------!', flush=True)
+    files = glob(inputs.folder+'/*.csv')
+    print('Total files: ', len(files), flush=True)
+
+    print('Reading files ...', flush=True)
+    _ = pool.map(map_files, files)
+    print('\tdone!', flush=True)
+
+    # 取出光谱数据
+    # data = [d[0] for d in _]
+    # for x in files:
+    #     map_files(files)
+
+    data = [d[0] for d in _]
+    data = xr.concat(data, dim='Measures')
+    # data.to_netcdf('C:\\fhz\\neibuwenjian\\tangchao\_Data\\temp.nc')
+
+    # 取出元数据
+    header = _[0][1]
+    _ = []
+    header = pd.DataFrame(header).set_index('Model')
+    # print(header, flush=True)
+
+    # sif算法处理
+    processing(inputs.standard_sif, inputs.folder, inputs.outfile,
+               inputs.pars, data, header, inputs.sky_p, inputs.method)
+
+
+if __name__ == '__main__':
+    main()
diff --git a/standard_sif.csv b/standard_sif.csv
new file mode 100644
index 0000000..76cb7dc
--- /dev/null
+++ b/standard_sif.csv
@@ -0,0 +1,1022 @@
+0,1
+649.223,0.00423238
+649.393,0.00441355
+649.563,0.00459897
+649.733,0.00478934
+649.903,0.00498536
+650.073,0.00518771
+650.243,0.00539664
+650.413,0.00561192
+650.583,0.00583325
+650.753,0.00606034
+650.923,0.00629292
+651.092,0.00652932
+651.262,0.00677241
+651.432,0.00702125
+651.602,0.00727628
+651.771,0.00753637
+651.941,0.00780503
+652.111,0.0080811
+652.28,0.00836216
+652.45,0.00865009
+652.619,0.00894007
+652.789,0.00923408
+652.958,0.00952719
+653.127,0.00981998
+653.297,0.0101154
+653.466,0.0104124
+653.635,0.0107151
+653.805,0.0110279
+653.974,0.0113494
+654.143,0.0116835
+654.312,0.0120289
+654.482,0.0123853
+654.651,0.0127459
+654.82,0.0131104
+654.989,0.0134762
+655.158,0.0138415
+655.327,0.0142071
+655.496,0.0145748
+655.665,0.0149462
+655.833,0.0153207
+656.002,0.0157046
+656.171,0.0160967
+656.34,0.0164963
+656.509,0.0169017
+656.677,0.0173092
+656.846,0.0177222
+657.015,0.0181368
+657.183,0.0185498
+657.352,0.0189668
+657.52,0.0193839
+657.689,0.0198073
+657.857,0.020233
+658.026,0.0206672
+658.194,0.0211056
+658.363,0.0215533
+658.531,0.0220049
+658.699,0.0224631
+658.867,0.0229277
+659.036,0.0234014
+659.204,0.0238784
+659.372,0.0243608
+659.54,0.0248479
+659.708,0.025339
+659.876,0.0258334
+660.044,0.0263304
+660.212,0.0268296
+660.38,0.0273309
+660.548,0.0278344
+660.716,0.0283401
+660.884,0.028848
+661.052,0.0293581
+661.22,0.0298707
+661.388,0.0303864
+661.555,0.0309026
+661.723,0.031426
+661.891,0.0319544
+662.058,0.0324849
+662.226,0.0330246
+662.393,0.0335673
+662.561,0.0341202
+662.728,0.0346769
+662.896,0.0352446
+663.063,0.0358168
+663.231,0.0363999
+663.398,0.0369846
+663.566,0.0375758
+663.733,0.038164
+663.9,0.0387505
+664.067,0.0393332
+664.235,0.0399154
+664.402,0.0404932
+664.569,0.041073
+664.736,0.0416578
+664.903,0.0422507
+665.07,0.0428547
+665.237,0.0434707
+665.404,0.0440976
+665.571,0.0447337
+665.738,0.0453776
+665.905,0.0460279
+666.072,0.0466829
+666.238,0.0473382
+666.405,0.0480015
+666.572,0.048669
+666.739,0.0493409
+666.905,0.0500132
+667.072,0.0506942
+667.239,0.0513802
+667.405,0.0520682
+667.572,0.0527676
+667.738,0.0534709
+667.905,0.0541877
+668.071,0.0549104
+668.237,0.0556432
+668.404,0.0563892
+668.57,0.0571379
+668.736,0.0578923
+668.903,0.0586554
+669.069,0.0594167
+669.235,0.0601809
+669.401,0.0609499
+669.567,0.0617261
+669.734,0.0625164
+669.9,0.0633134
+670.066,0.0641241
+670.232,0.0649489
+670.398,0.0657863
+670.564,0.0666344
+670.729,0.0674863
+670.895,0.0683505
+671.061,0.0692201
+671.227,0.070095
+671.393,0.0709772
+671.559,0.0718688
+671.724,0.0727666
+671.89,0.0736836
+672.056,0.0746167
+672.221,0.0755607
+672.387,0.0765259
+672.552,0.0774992
+672.718,0.0784908
+672.883,0.0794875
+673.049,0.0804997
+673.214,0.0815147
+673.379,0.0825395
+673.545,0.0835812
+673.71,0.0846283
+673.875,0.0856881
+674.041,0.0867681
+674.206,0.0878542
+674.371,0.0889483
+674.536,0.0900451
+674.701,0.0911395
+674.866,0.0922264
+675.031,0.0933006
+675.196,0.0943616
+675.361,0.0954194
+675.526,0.0964856
+675.691,0.0975718
+675.856,0.0986896
+676.021,0.0998506
+676.186,0.101065
+676.35,0.102329
+676.515,0.103665
+676.68,0.105068
+676.844,0.106534
+677.009,0.108087
+677.174,0.109719
+677.338,0.111416
+677.503,0.113192
+677.667,0.11502
+677.832,0.116916
+677.996,0.118852
+678.161,0.120843
+678.325,0.122852
+678.489,0.124875
+678.654,0.126909
+678.818,0.128914
+678.982,0.130886
+679.146,0.132814
+679.311,0.134712
+679.475,0.136566
+679.639,0.138399
+679.803,0.140218
+679.967,0.142034
+680.131,0.143855
+680.295,0.145681
+680.459,0.147508
+680.623,0.149333
+680.787,0.151152
+680.95,0.152951
+681.114,0.154748
+681.278,0.156534
+681.442,0.158315
+681.605,0.160085
+681.769,0.161869
+681.933,0.163663
+682.096,0.16546
+682.26,0.167278
+682.423,0.169083
+682.587,0.170887
+682.75,0.172656
+682.914,0.174402
+683.077,0.176092
+683.241,0.177747
+683.404,0.179355
+683.567,0.180938
+683.73,0.182507
+683.894,0.18408
+684.057,0.18565
+684.22,0.187229
+684.383,0.188806
+684.546,0.190374
+684.709,0.19192
+684.872,0.193437
+685.035,0.194913
+685.198,0.196343
+685.361,0.197729
+685.524,0.199074
+685.687,0.200382
+685.85,0.201655
+686.013,0.202898
+686.175,0.204104
+686.338,0.205283
+686.501,0.206423
+686.664,0.207516
+686.826,0.208552
+686.989,0.209536
+687.151,0.210452
+687.314,0.211312
+687.476,0.212107
+687.639,0.21285
+687.801,0.213533
+687.964,0.214167
+688.126,0.214746
+688.288,0.215278
+688.451,0.215766
+688.613,0.216207
+688.775,0.216608
+688.937,0.216969
+689.1,0.217296
+689.262,0.217585
+689.424,0.217836
+689.586,0.218048
+689.748,0.218219
+689.91,0.218347
+690.072,0.218432
+690.234,0.218473
+690.396,0.218471
+690.558,0.218427
+690.719,0.218342
+690.881,0.218217
+691.043,0.218053
+691.205,0.217852
+691.366,0.217624
+691.528,0.217374
+691.69,0.217111
+691.851,0.216845
+692.013,0.216579
+692.174,0.216325
+692.336,0.216085
+692.497,0.215869
+692.659,0.21568
+692.82,0.215525
+692.982,0.21541
+693.143,0.21534
+693.304,0.215313
+693.466,0.215325
+693.627,0.215373
+693.788,0.215453
+693.949,0.215561
+694.11,0.215693
+694.271,0.215845
+694.432,0.216015
+694.593,0.216199
+694.754,0.216392
+694.915,0.216593
+695.076,0.216796
+695.237,0.217001
+695.398,0.217208
+695.559,0.217416
+695.72,0.217627
+695.881,0.21784
+696.041,0.218055
+696.202,0.218277
+696.363,0.218511
+696.523,0.218763
+696.684,0.219043
+696.844,0.219356
+697.005,0.219712
+697.165,0.220111
+697.326,0.220554
+697.486,0.221029
+697.647,0.221533
+697.807,0.222054
+697.967,0.222589
+698.128,0.223133
+698.288,0.22368
+698.448,0.224233
+698.608,0.224792
+698.769,0.225364
+698.929,0.225941
+699.089,0.226528
+699.249,0.227124
+699.409,0.227728
+699.569,0.228338
+699.729,0.228952
+699.889,0.22957
+700.049,0.23019
+700.209,0.230812
+700.368,0.231437
+700.528,0.232079
+700.688,0.23274
+700.848,0.233424
+701.007,0.234132
+701.167,0.234877
+701.327,0.235654
+701.486,0.236454
+701.646,0.237285
+701.805,0.238133
+701.965,0.239006
+702.124,0.239892
+702.284,0.240803
+702.443,0.241731
+702.603,0.24269
+702.762,0.243671
+702.921,0.244684
+703.08,0.245731
+703.24,0.246819
+703.399,0.24793
+703.558,0.249064
+703.717,0.250216
+703.876,0.251383
+704.035,0.252559
+704.194,0.253741
+704.353,0.254927
+704.512,0.256119
+704.671,0.257314
+704.83,0.258513
+704.989,0.259717
+705.148,0.260924
+705.307,0.262135
+705.465,0.263345
+705.624,0.264569
+705.783,0.265802
+705.942,0.267044
+706.1,0.268289
+706.259,0.269555
+706.417,0.270829
+706.576,0.272128
+706.734,0.273439
+706.893,0.274782
+707.051,0.276142
+707.21,0.277538
+707.368,0.278952
+707.526,0.280396
+707.685,0.281877
+707.843,0.283378
+708.001,0.28491
+708.159,0.286471
+708.318,0.288074
+708.476,0.289697
+708.634,0.291351
+708.792,0.293035
+708.95,0.294751
+709.108,0.296497
+709.266,0.298274
+709.424,0.300081
+709.582,0.301916
+709.739,0.303766
+709.897,0.305655
+710.055,0.307569
+710.213,0.30951
+710.371,0.311481
+710.528,0.313472
+710.686,0.315515
+710.844,0.317598
+711.001,0.319714
+711.159,0.32189
+711.316,0.324094
+711.474,0.326346
+711.631,0.328611
+711.789,0.330911
+711.946,0.333208
+712.103,0.335511
+712.261,0.337832
+712.418,0.34014
+712.575,0.342448
+712.733,0.344773
+712.89,0.347082
+713.047,0.349391
+713.204,0.3517
+713.361,0.354011
+713.518,0.356328
+713.675,0.358652
+713.832,0.360987
+713.989,0.363335
+714.146,0.365701
+714.303,0.368094
+714.46,0.37053
+714.617,0.373019
+714.774,0.375576
+714.93,0.378195
+715.087,0.380924
+715.244,0.383744
+715.4,0.386624
+715.557,0.389589
+715.714,0.392607
+715.87,0.395646
+716.027,0.398732
+716.183,0.401814
+716.34,0.404918
+716.496,0.407992
+716.652,0.411042
+716.809,0.414076
+716.965,0.417042
+717.121,0.419949
+717.278,0.422819
+717.434,0.425625
+717.59,0.42839700000000003
+717.746,0.431146
+717.902,0.433883
+718.059,0.436635
+718.215,0.439379
+718.371,0.442146
+718.527,0.444949
+718.683,0.447801
+718.839,0.450714
+718.994,0.453683
+719.15,0.456756
+719.306,0.459909
+719.462,0.463134
+719.618,0.466424
+719.773,0.469748
+719.929,0.473142
+720.085,0.476575
+720.24,0.48002
+720.396,0.483515
+720.552,0.487031
+720.707,0.490541
+720.863,0.494084
+721.018,0.497609
+721.173,0.501135
+721.329,0.504683
+721.484,0.508205
+721.64,0.511746
+721.795,0.515261
+721.95,0.518769
+722.105,0.522272
+722.261,0.525789
+722.416,0.529274
+722.571,0.532748
+722.726,0.536211
+722.881,0.539661
+723.036,0.543097
+723.191,0.546517
+723.346,0.549925
+723.501,0.553321
+723.656,0.556708
+723.811,0.560087
+723.965,0.563439
+724.12,0.566809
+724.275,0.570177
+724.43,0.573548
+724.584,0.576901
+724.739,0.580282
+724.894,0.583674
+725.048,0.587056
+725.203,0.590474
+725.357,0.593885
+725.512,0.597334
+725.666,0.600777
+725.821,0.604258
+725.975,0.607734
+726.13,0.61125
+726.284,0.614765
+726.438,0.618304
+726.592,0.621873
+726.747,0.625498
+726.901,0.629138
+727.055,0.63282
+727.209,0.636544
+727.363,0.640303
+727.517,0.644093
+727.671,0.647905
+727.825,0.651735
+727.979,0.655576
+728.133,0.659421
+728.287,0.663268
+728.441,0.667112
+728.595,0.670951
+728.749,0.674782
+728.902,0.678577
+729.056,0.682382
+729.21,0.686172
+729.363,0.689922
+729.517,0.693683
+729.671,0.697431
+729.824,0.701143
+729.978,0.704869
+730.131,0.70856
+730.285,0.712269
+730.438,0.715951
+730.592,0.719656
+730.745,0.723341
+730.898,0.727033
+731.051,0.730735
+731.205,0.734473
+731.358,0.738193
+731.511,0.741917
+731.664,0.745641
+731.817,0.749361
+731.971,0.753098
+732.124,0.756798
+732.277,0.760479
+732.43,0.764135
+732.583,0.767757
+732.736,0.771338
+732.888,0.774848
+733.041,0.778326
+733.194,0.781743
+733.347,0.785102
+733.5,0.788402
+733.652,0.791624
+733.805,0.794813
+733.958,0.797949
+734.11,0.801012
+734.263,0.80405
+734.416,0.80705
+734.568,0.810002
+734.721,0.812951
+734.873,0.815867
+735.025,0.818779
+735.178,0.821709
+735.33,0.824617
+735.482,0.827518
+735.635,0.830427
+735.787,0.833303
+735.939,0.83616
+736.091,0.838994
+736.244,0.841824
+736.396,0.844613
+736.548,0.847382
+736.7,0.85013
+736.852,0.852859
+737.004,0.855571
+737.156,0.858265
+737.308,0.860935
+737.46,0.863574
+737.612,0.866176
+737.763,0.868717
+737.915,0.871223
+738.067,0.873671
+738.219,0.876053
+738.37,0.878345
+738.522,0.880567
+738.674,0.882695
+738.825,0.884708
+738.977,0.886621
+739.128,0.888403
+739.28,0.890078
+739.431,0.891633
+739.583,0.893091
+739.734,0.894442
+739.885,0.8957
+740.037,0.896879
+740.188,0.897971
+740.339,0.89899
+740.491,0.899951
+740.642,0.900849
+740.793,0.901699
+740.944,0.902509
+741.095,0.903285
+741.246,0.904027
+741.397,0.904729
+741.548,0.905384
+741.699,0.905985
+741.85,0.906527
+742.001,0.907003
+742.152,0.907406
+742.303,0.907732
+742.453,0.907975
+742.604,0.908132
+742.755,0.908199
+742.906,0.908169
+743.056,0.90804
+743.207,0.907808
+743.357,0.907478
+743.508,0.907049
+743.659,0.906524
+743.809,0.905911
+743.96,0.905202
+744.11,0.904412
+744.26,0.903538
+744.411,0.902581
+744.561,0.901558
+744.711,0.900468
+744.862,0.899309
+745.012,0.898101
+745.162,0.896841
+745.312,0.895526
+745.462,0.894153
+745.612,0.892719
+745.762,0.891219
+745.912,0.889652
+746.062,0.888014
+746.212,0.886301
+746.362,0.884513
+746.512,0.882645
+746.662,0.880697
+746.812,0.878666
+746.962,0.87655
+747.111,0.874362
+747.261,0.872075
+747.411,0.869706
+747.56,0.867278
+747.71,0.86476
+747.86,0.862173
+748.009,0.859539
+748.159,0.856826
+748.308,0.854073
+748.458,0.851245
+748.607,0.848384
+748.757,0.845452
+748.906,0.842491
+749.055,0.839485
+749.205,0.836417
+749.354,0.83333
+749.503,0.830209
+749.652,0.827057
+749.801,0.823879
+749.95,0.820678
+750.1,0.817437
+750.249,0.8142
+750.398,0.810945
+750.547,0.807671
+750.696,0.804378
+750.845,0.801066
+750.993,0.797757
+751.142,0.794406
+751.291,0.791036
+751.44,0.787649
+751.589,0.784247
+751.737,0.780855
+751.886,0.777428
+752.035,0.773992
+752.183,0.77057
+752.332,0.767116
+752.48,0.763675
+752.629,0.760201
+752.777,0.756739
+752.926,0.753241
+753.074,0.749756
+753.223,0.746236
+753.371,0.742735
+753.519,0.739232
+753.668,0.73571
+753.816,0.73222
+753.964,0.728743
+754.112,0.725284
+754.261,0.721823
+754.409,0.71841
+754.557,0.715026
+754.705,0.711674
+754.853,0.708357
+755.001,0.705078
+755.149,0.70184
+755.297,0.698634
+755.445,0.695454
+755.592,0.692311
+755.74,0.689154
+755.888,0.685996
+756.036,0.682828
+756.184,0.679642
+756.331,0.676451
+756.479,0.673206
+756.627,0.669919
+756.774,0.666605
+756.922,0.663213
+757.069,0.65978
+757.217,0.656258
+757.364,0.652694
+757.512,0.649045
+757.659,0.645361
+757.806,0.64162
+757.954,0.637798
+758.101,0.633951
+758.248,0.630055
+758.395,0.626114
+758.543,0.622105
+758.69,0.618085
+758.837,0.614031
+758.984,0.609946
+759.131,0.605836
+759.278,0.601708
+759.425,0.597572
+759.572,0.593438
+759.719,0.589315
+759.866,0.585213
+760.013,0.581141
+760.159,0.577138
+760.306,0.573158
+760.453,0.56924
+760.6,0.565397
+760.746,0.561663
+760.893,0.558
+761.04,0.554443
+761.186,0.551021
+761.333,0.547677
+761.479,0.544446
+761.626,0.541275
+761.772,0.538194
+761.919,0.535154
+762.065,0.532183
+762.211,0.529253
+762.358,0.526336
+762.504,0.523464
+762.65,0.520611
+762.797,0.517749
+762.943,0.514909
+763.089,0.512066
+763.235,0.509216
+763.381,0.50636
+763.527,0.503498
+763.673,0.500632
+763.819,0.497762
+763.965,0.494889
+764.111,0.492014
+764.257,0.489138
+764.403,0.486265
+764.549,0.483397
+764.694,0.480557
+764.84,0.477707
+764.986,0.474871
+765.132,0.47205
+765.277,0.469261
+765.423,0.466462
+765.568,0.463685
+765.714,0.46089
+765.859,0.458111
+766.005,0.455304
+766.15,0.452504
+766.296,0.449665
+766.441,0.446821
+766.587,0.443925
+766.732,0.441011
+766.877,0.438051
+767.022,0.43504
+767.168,0.431953
+767.313,0.428842
+767.458,0.425696
+767.603,0.422525
+767.748,0.41934
+767.893,0.41615
+768.038,0.412967
+768.183,0.409798
+768.328,0.406651
+768.473,0.40353
+768.618,0.400441
+768.763,0.39739
+768.908,0.394383
+769.052,0.391445
+769.197,0.388538
+769.342,0.385679
+769.487,0.382864
+769.631,0.380107
+769.776,0.377369
+769.92,0.374681
+770.065,0.372003
+770.21,0.369348
+770.354,0.366725
+770.499,0.364089
+770.643,0.361468
+770.787,0.358835
+770.932,0.356163
+771.076,0.35348
+771.22,0.350768
+771.365,0.348012
+771.509,0.345257
+771.653,0.342489
+771.797,0.339713
+771.941,0.336937
+772.085,0.334164
+772.229,0.331401
+772.373,0.328649
+772.517,0.325914
+772.661,0.323198
+772.805,0.320504
+772.949,0.317838
+773.093,0.3152
+773.237,0.312596
+773.381,0.310025
+773.524,0.307507
+773.668,0.305008
+773.812,0.302549
+773.955,0.300147
+774.099,0.297771
+774.243,0.295437
+774.386,0.293163
+774.53,0.290918
+774.673,0.288731
+774.817,0.286575
+774.96,0.284478
+775.103,0.282427
+775.247,0.280403
+775.39,0.27843
+775.533,0.276488
+775.677,0.274559
+775.82,0.272664
+775.963,0.270784
+776.106,0.268915
+776.249,0.267051
+776.392,0.265186
+776.535,0.263313
+776.678,0.261427
+776.821,0.259521
+776.964,0.257591
+777.107,0.25563
+777.25,0.253642
+777.393,0.251636
+777.536,0.24962
+777.679,0.2476
+777.821,0.245601
+777.964,0.2436
+778.107,0.241622
+778.249,0.239685
+778.392,0.237771
+778.535,0.2359
+778.677,0.234092
+778.82,0.232328
+778.962,0.23064
+779.105,0.229011
+779.247,0.227461
+779.389,0.225971
+779.532,0.224527
+779.674,0.22314
+779.816,0.221794
+779.959,0.220473
+780.101,0.21919
+780.243,0.217927
+780.385,0.21668
+780.527,0.215442
+780.669,0.214206
+780.811,0.212967
+780.953,0.211717
+781.095,0.210451
+781.237,0.209169
+781.379,0.207872
+781.521,0.206564
+781.663,0.205247
+781.805,0.203924
+781.947,0.202596
+782.088,0.201276
+782.23,0.199949
+782.372,0.198626
+782.513,0.197319
+782.655,0.196013
+782.797,0.19472
+782.938,0.193452
+783.08,0.192193
+783.221,0.190963
+783.363,0.189744
+783.504,0.188552
+783.645,0.187378
+783.787,0.186214
+783.928,0.185075
+784.069,0.183953
+784.211,0.182839
+784.352,0.181749
+784.493,0.180672
+784.634,0.179608
+784.775,0.178557
+784.916,0.177515
+785.057,0.176484
+785.198,0.175463
+785.339,0.174452
+785.48,0.173454
+785.621,0.172469
+785.762,0.1715
+785.903,0.170546
+786.044,0.16961
+786.185,0.168691
+786.325,0.167794
+786.466,0.166903
+786.607,0.166022
+786.747,0.165155
+786.888,0.164287
+787.029,0.163422
+787.169,0.162564
+787.31,0.161699
+787.45,0.160837
+787.591,0.159966
+787.731,0.159095
+787.871,0.158217
+788.012,0.157324
+788.152,0.156426
+788.292,0.155516
+788.433,0.154583
+788.573,0.153637
+788.713,0.152668
+788.853,0.151674
+788.993,0.150652
+789.133,0.1496
+789.273,0.148523
+789.413,0.147428
+789.553,0.146322
+789.693,0.145212
+789.833,0.144106
+789.973,0.14301
+790.113,0.141931
+790.253,0.140877
+790.393,0.139856
+790.532,0.138883
+790.672,0.137951
+790.812,0.137074
+790.951,0.136268
+791.091,0.135526
+791.231,0.134855
+791.37,0.13425
+791.51,0.133697
+791.649,0.133196
+791.789,0.132733
+791.928,0.132309
+792.067,0.131912
+792.207,0.131533
+792.346,0.131173
+792.485,0.130821
+792.625,0.130471
+792.764,0.13012
+792.903,0.12976
+793.042,0.129385
+793.181,0.128992
+793.32,0.128583
+793.459,0.128158
+793.598,0.127721
+793.737,0.127273
+793.876,0.126815
+794.015,0.126349
+794.154,0.125878
+794.293,0.125404
+794.432,0.124927
+794.571,0.12445
+794.709,0.123979
+794.848,0.123508
+794.987,0.123043
+795.125,0.122588
+795.264,0.122137
+795.403,0.121691
+795.541,0.121252
+795.68,0.120811
+795.818,0.120376
+795.957,0.119936
+796.095,0.119498
+796.233,0.119059
+796.372,0.118614
+796.51,0.118172
+796.648,0.117729
+796.787,0.117283
+796.925,0.11684
+797.063,0.116398
+797.201,0.115957
+797.339,0.115516
+797.477,0.115075
+797.615,0.114634
+797.753,0.114192
+797.891,0.11375
+798.029,0.113307
+798.167,0.112862
+798.305,0.112417
+798.443,0.111973
+798.581,0.11153
+798.719,0.111089
+798.856,0.110654
+798.994,0.110219
+799.132,0.109788
+799.269,0.109364
+799.407,0.108938
+799.545,0.108513
+799.682,0.108091
+799.82,0.107663
+799.957,0.107235
+800.095,0.1068
+800.232,0.106363
+800.369,0.105924
+800.507,0.10548
+800.644,0.105038
+800.781,0.104598
+800.919,0.104157
+801.056,0.103723
+801.193,0.103293
+801.33,0.102866
+801.467,0.102441
+801.604,0.102019
+801.741,0.101597
+801.878,0.101176
+802.015,0.100754
+802.152,0.100331
+802.289,0.0999067
+802.426,0.0994819
+802.563,0.0990567
+802.7,0.0986313
+802.837,0.0982058
+802.973,0.0977837
+803.11,0.097359
+803.247,0.0969347
+803.383,0.0965136
+803.52,0.0960894
+803.656,0.0956681
+803.793,0.0952432
+803.93,0.0948177
+804.066,0.0943945
+804.202,0.0939707
+804.339,0.0935436
+804.475,0.0931202
+804.612,0.0926946
+804.749,0.0922707
+804.886,0.0918489
+805.023,0.0914299
+805.16,0.0910138
+805.297,0.0906004
+805.434,0.090189