[pktools] 107/375: removed vito related code from master
Bas Couwenberg
sebastic at xs4all.nl
Wed Dec 3 21:54:04 UTC 2014
This is an automated email from the git hooks/post-receive script.
sebastic-guest pushed a commit to branch upstream-master
in repository pktools.
commit 7872221aad6bbb7cac8488e8649668e61bca225f
Author: user <user at osgeolive.(none)>
Date: Tue May 14 21:22:23 2013 +0200
removed vito related code from master
---
Makefile.am | 1 -
configure.ac | 2 -
src/apps/Makefile.am | 7 +-
src/apps/pkgetchandelier.cc | 361 -------------
src/apps/pkgetchandelier.h | 7 -
src/apps/pksensormodel.cc | 668 -------------------------
src/apps/pksetchandelier.cc | 224 ---------
src/base/Makefile.am | 5 +-
src/base/PointData.cc | 34 --
src/base/PointData.h | 51 --
src/models/Makefile.am | 40 --
src/models/Prospect.cc | 31 --
src/models/Prospect.h | 58 ---
src/models/SensorModel.h | 399 ---------------
src/models/dataSpec_P5B.f90 | 1139 ------------------------------------------
src/models/pktestProspect.cc | 76 ---
src/models/prospect_5B.f90 | 132 -----
src/models/tav_abs.f90 | 60 ---
18 files changed, 3 insertions(+), 3292 deletions(-)
diff --git a/Makefile.am b/Makefile.am
index 0fc0a98..d3bf4bd 100755
--- a/Makefile.am
+++ b/Makefile.am
@@ -5,6 +5,5 @@ ACLOCAL_AMFLAGS = -I m4
SUBDIRS = src/base \
src/imageclasses \
src/algorithms \
- src/models \
$(FILECLASSES_OPT) \
src/apps
diff --git a/configure.ac b/configure.ac
index 0a045a0..4a575a9 100644
--- a/configure.ac
+++ b/configure.ac
@@ -17,7 +17,6 @@ AC_MSG_WARN([$GDAL_OGR_ENABLED])
AC_PROG_CXX
AC_PROG_FC
LT_INIT
-#AC_PROG_RANLIB
# check if the source folder is correct
AC_CONFIG_SRCDIR([src/apps/pkinfo.cc])
@@ -98,7 +97,6 @@ AC_CONFIG_FILES([
Makefile
src/base/Makefile
src/algorithms/Makefile
-src/models/Makefile
src/imageclasses/Makefile
src/fileclasses/Makefile
src/apps/Makefile
diff --git a/src/apps/Makefile.am b/src/apps/Makefile.am
index f717c7e..7dd7113 100644
--- a/src/apps/Makefile.am
+++ b/src/apps/Makefile.am
@@ -28,10 +28,7 @@ pklas2img_SOURCES = pklas2img.cc
pklas2img_LDADD = -llas $(AM_LDFLAGS)
endif
if USE_NLOPT
-bin_PROGRAMS += pksensormodel pkopt_svm
-#pksensormodel_LDFLAGS = -O1
-pksensormodel_SOURCES = $(top_srcdir)/src/models/SensorModel.h pksensormodel.h pksensormodel.cc
-pksensormodel_LDADD = $(GSL_LIBS) $(AM_LDFLAGS) -lnlopt -lm -larmadillo
+bin_PROGRAMS += pkopt_svm
pkopt_svm_SOURCES = $(top_srcdir)/src/algorithms/OptFactory.h pkclassify_nn.h pkopt_svm.cc
pkopt_svm_LDADD = $(GSL_LIBS) $(AM_LDFLAGS) -lnlopt
endif
@@ -67,6 +64,4 @@ pkascii2img_SOURCES = pkascii2img.cc
pkascii2ogr_SOURCES = pkascii2ogr.cc
pkeditogr_SOURCES = pkeditogr.cc
pkeditogr_LDADD = $(AM_LDFLAGS)
-#pkxcorimg_SOURCES= pkxcorimg.cc
-#pkgeom_SOURCES = pkgeom.cc
###############################################################################
diff --git a/src/apps/pkgetchandelier.cc b/src/apps/pkgetchandelier.cc
deleted file mode 100644
index 46b4a15..0000000
--- a/src/apps/pkgetchandelier.cc
+++ /dev/null
@@ -1,361 +0,0 @@
-/**********************************************************************
-pkgetchandelier.cc: program to optimize model parameters for brdf correction
-Copyright (C) 2008-2012 Pieter Kempeneers
-
-This file is part of pktools
-
-pktools is free software: you can redistribute it and/or modify
-it under the terms of the GNU General Public License as published by
-the Free Software Foundation, either version 3 of the License, or
-(at your option) any later version.
-
-pktools is distributed in the hope that it will be useful,
-but WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-GNU General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with pktools. If not, see <http://www.gnu.org/licenses/>.
-***********************************************************************/
-#include <assert.h>
-#include <string>
-#include <iostream>
-#include <nlopt.hpp>
-#include "base/PointData.h"
-#include "algorithms/StatFactory.h"
-#include "imageclasses/ImgReaderOgr.h"
-#include "Optionpk.h"
-#include "pkgetchandelier.h"
-
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
-
-double objFunction(const std::vector<double> &x, std::vector<double> &grad, void *my_func_data){
- assert(grad.empty());
- double rmse=0;
- vector<vector<PointData> > *blockData=reinterpret_cast<vector< vector <PointData> > *> (my_func_data);
- for(unsigned int ipoint=0;ipoint<blockData->size();++ipoint){
- for(unsigned int ifile1=0;ifile1<blockData->at(ipoint).size()-1;++ifile1){
- double reflectance1=0;
- if(blockData->at(ipoint)[ifile1].getImage()<0)
- reflectance1=blockData->at(ipoint)[ifile1].getReflectance();
- else{
- State state1;
- state1.k=x[blockData->at(ipoint)[ifile1].getImage()*4];
- state1.e=x[blockData->at(ipoint)[ifile1].getImage()*4+1];
- state1.a=x[blockData->at(ipoint)[ifile1].getImage()*4+2];
- state1.haze=x[blockData->at(ipoint)[ifile1].getImage()*4+3];
- reflectance1=blockData->at(ipoint)[ifile1].correctReflectance(state1);
- }
- for(unsigned int ifile2=ifile1;ifile2<blockData->at(ipoint).size();++ifile2){
- double reflectance2=0;
- if(blockData->at(ipoint)[ifile2].getImage()<0)
- reflectance2=blockData->at(ipoint)[ifile2].getReflectance();
- else{
- State state2;
- state2.k=x[blockData->at(ipoint)[ifile2].getImage()*4];
- state2.e=x[blockData->at(ipoint)[ifile2].getImage()*4+1];
- state2.a=x[blockData->at(ipoint)[ifile2].getImage()*4+2];
- state2.haze=x[blockData->at(ipoint)[ifile2].getImage()*4+3];
- reflectance2=blockData->at(ipoint)[ifile2].correctReflectance(state2);
- }
- // if(reflectance2>2*reflectance1||reflectance2<reflectance1/2)
- // continue;
- double diff=(reflectance1-reflectance2)*(reflectance1-reflectance2);
- if(PointData::m_residual>0){
- if(sqrt(diff/reflectance1/reflectance2)<PointData::m_residual)
- rmse+=diff;
- }
- else
- rmse+=diff;
- }
- }
- }
- rmse=sqrt(rmse/blockData->size());
- // std::cout << "difference: " << diff << std::endl;
- return rmse;
-}
-
-int main(int argc, char *argv[])
-{
- std::string versionString="version ";
- versionString+=VERSION;
- versionString+=", Copyright (C) 2008-2012 Pieter Kempeneers.\n\
- This program comes with ABSOLUTELY NO WARRANTY; for details type use option -h.\n\
- This is free software, and you are welcome to redistribute it\n\
- under certain conditions; use option --license for details.";
- Optionpk<bool> version_opt("\0","version",versionString,false);
- Optionpk<bool> license_opt("lic","license","show license information",false);
- Optionpk<bool> help_opt("h","help","shows this help info",false);
- Optionpk<bool> todo_opt("\0","todo","",false);
- Optionpk<string> input_opt("i", "input", "Reflectance imageinput (ogr vector) files","");
- // Optionpk<string> sample_opt("s", "sample", "sample (ogr vector) file","");
- Optionpk<double> grid_opt("grid", "grid", "grid for tiepoints (in m). Smaller values provide more tiepoints",100);
- Optionpk<string> geom_opt("g", "geom", "Geometry image files","");
- Optionpk<unsigned short> band_opt("b", "band", "field name of Reflectance band in reflectance file",0);
- Optionpk<unsigned short> sza_opt("sza", "sza", "band number (starting from 0) for Sun Zenith Angle in geometry image file",3);
- // Optionpk<unsigned short> vza_opt("vza", "vza", "band number (starting from 0) for Sun Zenith Angle in geometry image file",5);
- Optionpk<unsigned short> saa_opt("saa", "saa", "band number (starting from 0) for Sun Azimuth Angle in geometry image file",6);
- Optionpk<unsigned short> vaa_opt("vaa", "vaa", "band number (starting from 0) for View Azimuth Angle in geometry image file",7);
- Optionpk<double> tau_opt("tau", "tau", "Optical depth",0.2);
- Optionpk<double> deltaT_opt("deltaT", "deltaT", "Exposure time",1.0);
- Optionpk<double> residual_opt("res", "residual", "residual error for filtering outlier tiepoints: maximum relative error to take tie point into account (use 0 if not filtering is required)",0.0);
- Optionpk<unsigned int> maxit_opt("maxit","maxit","maximum number of iterations",500);
- Optionpk<double> tolerance_opt("tol","tolerance","relative tolerance for stopping criterion",0.0001);
- Optionpk<double> init_opt("is","init","initial state vector: -is k -is e -is a -s haze",0);
- Optionpk<double> lb_opt("lb","lb","lower bounds for k,e and a: use -lb k -lb e -lb a -lb haze",0);
- Optionpk<double> ub_opt("ub","ub","upper bounds for k,e and a: use -ub k -ub e -ub a -ub haze",0);
- Optionpk<short> dim_opt("dim","dim","window size to calculate mean reflectance (use odd number, e.g., 3, 5, 7)",3);
- Optionpk<double> var_opt("var","var","maximum variance in window to take tiepoint into account (use 0 if don't care)",0);
- Optionpk<string> mask_opt("m", "mask", "Mask image(s) (single mask for all input images or one mask for each input image", "");
- Optionpk<double> invalid_opt("t", "invalid", "Mask value where image is invalid.", 0);
- // Optionpk<bool> single_opt("\0","single","create virtual tiepoints for single images with no overlap",false);
- Optionpk<unsigned long int> mintp_opt("min", "min", "Minimum number of tiepoints", 1000);
- Optionpk<short> verbose_opt("v", "verbose", "verbose", 0);
-
- version_opt.retrieveOption(argc,argv);
- license_opt.retrieveOption(argc,argv);
- help_opt.retrieveOption(argc,argv);
- todo_opt.retrieveOption(argc,argv);
-
- input_opt.retrieveOption(argc,argv);
- // sample_opt.retrieveOption(argc,argv);
- grid_opt.retrieveOption(argc,argv);
- geom_opt.retrieveOption(argc,argv);
- band_opt.retrieveOption(argc,argv);
- sza_opt.retrieveOption(argc,argv);
- saa_opt.retrieveOption(argc,argv);
- vaa_opt.retrieveOption(argc,argv);
- tau_opt.retrieveOption(argc,argv);
- deltaT_opt.retrieveOption(argc,argv);
- residual_opt.retrieveOption(argc,argv);
- maxit_opt.retrieveOption(argc,argv);
- tolerance_opt.retrieveOption(argc,argv);
- init_opt.retrieveOption(argc,argv);
- lb_opt.retrieveOption(argc,argv);
- ub_opt.retrieveOption(argc,argv);
- dim_opt.retrieveOption(argc,argv);
- var_opt.retrieveOption(argc,argv);
- mask_opt.retrieveOption(argc,argv);
- invalid_opt.retrieveOption(argc,argv);
- // single_opt.retrieveOption(argc,argv);
- verbose_opt.retrieveOption(argc,argv);
- mintp_opt.retrieveOption(argc,argv);
-
- if(version_opt[0]||todo_opt[0]){
- cout << version_opt.getHelp() << endl;
- cout << "todo: " << todo_opt.getHelp() << endl;
- exit(0);
- }
- if(license_opt[0]){
- cout << Optionpk<bool>::getGPLv3License() << endl;
- exit(0);
- }
- if(help_opt[0]){
- cout << "usage: pktrirad -s tiepoints -i input [-i input] -g geom [-g geom] [OPTIONS]" << endl;
- exit(0);
- }
-
- assert(geom_opt.size()==input_opt.size());
-
- assert(init_opt.size()==4*input_opt.size());//k,e,a,haze
- assert(lb_opt.size()==4*input_opt.size());
- assert(ub_opt.size()==4*input_opt.size());
-
- PointData::m_tau=tau_opt[0];
- PointData::m_deltaT=deltaT_opt[0];
- PointData::m_residual=residual_opt[0];
- //construct list of points defined in sample file
- vector<vector< PointData> > blockData;//[point][image]
- vector<PointData> pdVector;//list of points (contains all image data for this tiepoint)
- //get union bounding box
- double maxLRX=0;
- double maxULY=0;
- double minULX=0;
- double minLRY=0;
- for(int ifile=0;ifile<input_opt.size();++ifile){
- double theULX, theULY, theLRX, theLRY;
- if(verbose_opt[0]>1)
- std::cout << "opening input file " << input_opt[ifile] << std::endl;
- ImgReaderGdal inputReader(input_opt[ifile]);
- inputReader.getBoundingBox(theULX,theULY,theLRX,theLRY);
- if(verbose_opt[0])
- std::cout << setprecision(12) << "--ulx=" << theULX << " --uly=" << theULY << " --lrx=" << theLRX << " --lry=" << theLRY << " ";
- if(!ifile){
- maxLRX=theLRX;
- maxULY=theULY;
- minULX=theULX;
- minLRY=theLRY;
- }
- else{
- maxLRX=(theLRX>maxLRX)?theLRX:maxLRX;
- maxULY=(theULY>maxULY)?theULY:maxULY;
- minULX=(theULX<minULX)?theULX:minULX;
- minLRY=(theLRY<minLRY)?theLRY:minLRY;
- }
- inputReader.close();
- }
- if(verbose_opt[0])
- std::cout << "union bounding box: " << setprecision(12) << "--ulx=" << minULX << " --uly=" << maxULY << " --lrx=" << maxLRX << " --lry=" << minLRY << std::endl;
-
- if(verbose_opt[0])
- cout << "number of mask images: " << mask_opt.size() << endl;
- vector<double> oldRowMask(mask_opt.size());
- int imask=0;
- ImgReaderGdal maskReader;
- if(mask_opt[0]!=""){
- if(mask_opt.size()!=input_opt.size())//single mask for all input images: open now
- maskReader.open(mask_opt[imask]);
- }
- for(double y=maxULY-grid_opt[0]/2.0;y>minLRY;y-=grid_opt[0]){
- for(double x=minULX+grid_opt[0]/2.0;x<maxLRX;x+=grid_opt[0]){
- for(int ifile=0;ifile<input_opt.size();++ifile){
- PointData pd;
- pd.setImage(ifile);
- if(verbose_opt[0]>1)
- std::cout << "opening input file " << input_opt[ifile] << std::endl;
- ImgReaderGdal inputReader(input_opt[ifile]);
- if(verbose_opt[0]>1)
- std::cout << "opening geom input file " << geom_opt[ifile] << std::endl;
- ImgReaderGdal geoReader(geom_opt[ifile]);
- double reflectance;
- double col,row;
- inputReader.geo2image(x,y,col,row);
- if(col<0||row<0||col>inputReader.nrOfCol()||row>inputReader.nrOfRow()){
- inputReader.close();
- geoReader.close();
- continue;
- }
- //begin mask
- if(mask_opt[0]!=""){
- bool masked=false;
- if(mask_opt.size()==input_opt.size())//one mask for each input image: open now
- maskReader.open(mask_opt[ifile]);
- double colMask,rowMask;
- maskReader.geo2image(x,y,colMask,rowMask);
- colMask=static_cast<int>(colMask);
- rowMask=static_cast<int>(rowMask);
- double maskValue=0;
- if(rowMask>=0&&rowMask<maskReader.nrOfRow()&&colMask>=0&&colMask<maskReader.nrOfCol()){
- try{
- maskReader.readData(maskValue,GDT_Float64,static_cast<int>(colMask),static_cast<int>(rowMask));
- }
- catch(string errorstring){
- cerr << errorstring << endl;
- exit(1);
- }
- if(maskValue==invalid_opt[0])
- masked=true;
- }
- if(mask_opt.size()==input_opt.size())//one mask for each input image: close now
- maskReader.close();
- if(masked){
- inputReader.close();
- geoReader.close();
- continue;
- }
- }
- //end mask
- statfactory::StatFactory stat;
- vector<double> windowBuffer;
- for(int windowJ=-dim_opt[0]/2;windowJ<(dim_opt[0]+1)/2;++windowJ){
- double j=row+windowJ;
- if(j<0||j>=inputReader.nrOfRow())
- continue;
- for(int windowI=-dim_opt[0]/2;windowI<(dim_opt[0]+1)/2;++windowI){
- double i=col+windowI;
- if(i<0||i>=inputReader.nrOfCol())
- continue;
- inputReader.readData(reflectance,GDT_Float64,i,j,band_opt[0]);
- windowBuffer.push_back(reflectance);
- }
- }
- double variance=0;
- if(windowBuffer.size()>1)
- stat.meanVar(windowBuffer,reflectance,variance);
- if(var_opt[0]>0){
- if(verbose_opt[0]>1)
- std::cout << "reflectance (" << windowBuffer.size() << " at " << col << "," << row << ") mean, var: " << reflectance << ", " << variance << std::endl;
- if(variance>var_opt[0]){
- inputReader.close();
- geoReader.close();
- continue;
- }
- }
- pd.setReflectance(reflectance);
- double sza;
- double saa;
- double vaa;
- geoReader.geo2image(x,y,col,row);
- geoReader.readData(sza,GDT_Float64,col,row,sza_opt[0]);
- geoReader.readData(saa,GDT_Float64,col,row,saa_opt[0]);
- geoReader.readData(vaa,GDT_Float64,col,row,vaa_opt[0]);
- double theRelativeAzimuth=saa-vaa;
- pd.setCosFi(theRelativeAzimuth);
- pd.setCosSZA(sza);
- double centreX,centreY;
- inputReader.getCentrePos(centreX,centreY);
- double ulx=inputReader.getUlx();
- double uly=inputReader.getUly();
- double r0=((ulx-centreX)*(ulx-centreX)+(uly-centreY)*(uly-centreY));
- double r=((x-centreX)*(x-centreX)+(y-centreY)*(y-centreY));
- pd.setR(r/r0);
- pdVector.push_back(pd);
- inputReader.close();
- geoReader.close();
- }
- if(pdVector.size()>1)//tiepoints must cover at least two images
- blockData.push_back(pdVector);
- // else if(single_opt[0]){
- // PointData pd=pdVector.back();
- // pd.setImage(-1);
- // pdVector.push_back(pd);
- // blockData.push_back(pdVector);
- // }
- pdVector.clear();
- if(mask_opt.size()!=input_opt.size())//single mask for all input images: close now
- maskReader.close();
- }
- }
-
- if(verbose_opt[0])
- std::cout << "Number of tiepoints: " << blockData.size() << std::endl;
- if(blockData.size()<mintp_opt[0]){
- std::cerr << "Error: Number of tiepoints=" << blockData.size() << " is smaller than " << mintp_opt[0] << std::endl;
- exit;
- }
-
- //todo: make nlopt::LN_COBYLA as a command line option
- //nlopt::opt opt(nlopt::LN_COBYLA,4*input_opt.size());//k,e,a,haze
- nlopt::opt optimizer(nlopt::LN_SBPLX,4*input_opt.size());//k,e,a,haze
-
- optimizer.set_min_objective(objFunction, &blockData);
-
- if(verbose_opt[0])
- std::cout << "set lower and upper bounds" << std::endl;
- optimizer.set_lower_bounds(lb_opt);
- optimizer.set_upper_bounds(ub_opt);
-
- if(verbose_opt[0])
- std::cout << "set stopping criteria" << std::endl;
- //set stopping criteria
- if(maxit_opt[0])
- optimizer.set_maxeval(maxit_opt[0]);
- else
- optimizer.set_xtol_rel(tolerance_opt[0]);
- double minf=0;
- std::vector<double> x=init_opt;
-
- if(verbose_opt[0])
- std::cout << "optimizing with " << optimizer.get_algorithm_name() << "..." << std::endl;
- optimizer.optimize(x, minf);
- for(int index=0;index<x.size();++index){
- std::cout << " -s " << x[index];
- if(index%4==3)
- std::cout << std::endl;
- }
- if(verbose_opt[0])
- std::cout << "minf: " << minf << std::endl;
-}
-
diff --git a/src/apps/pkgetchandelier.h b/src/apps/pkgetchandelier.h
deleted file mode 100644
index 889fcf6..0000000
--- a/src/apps/pkgetchandelier.h
+++ /dev/null
@@ -1,7 +0,0 @@
-#ifndef _PKGETCHANDELIER_
-#define _PKGETCHANDELIER_
-#include <math.h>
-#include <vector>
-
-double objFunction(const std::vector<double> &x, std::vector<double> &grad, void *my_func_data);
-#endif //_PKGETCHANDELIER_
diff --git a/src/apps/pksensormodel.cc b/src/apps/pksensormodel.cc
deleted file mode 100644
index 24594c7..0000000
--- a/src/apps/pksensormodel.cc
+++ /dev/null
@@ -1,668 +0,0 @@
-/**********************************************************************
-pksensormodel.cc: program to calculate geometric position based on row (sensor), col (sensor), roll, pitch, yaw and lens coordinates
-Copyright (C) 2008-2012 Pieter Kempeneers
-
-This file is part of pktools
-
-pktools is free software: you can redistribute it and/or modify
-it under the terms of the GNU General Public License as published by
-the Free Software Foundation, either version 3 of the License, or
-(at your option) any later version.
-
-pktools is distributed in the hope that it will be useful,
-but WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-GNU General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with pktools. If not, see <http://www.gnu.org/licenses/>.
-***********************************************************************/
-#include <iostream>
-#include <sstream>
-#include <fstream>
-#include <vector>
-#include <math.h>
-#include <nlopt.hpp>
-#include "base/Optionpk.h"
-#include "algorithms/OptFactory.h"
-#include "algorithms/StatFactory.h"
-#include "pksensormodel.h"
-
-double objFunction(const std::vector<double> &x, std::vector<double> &grad, void *my_func_data){
- assert(grad.empty());
- double error=0;
- DataModel *dm=reinterpret_cast<DataModel*> (my_func_data);
- arma::vec bc_att(3);
- bc_att(0)=x[0];
- bc_att(1)=x[1];
- bc_att(2)=x[2];
- dm->setBoresightAtt(bc_att);
- for(unsigned int ipoint=0;ipoint<dm->getSize();++ipoint){
- double e=dm->getDistGeo(ipoint);
- error+=e;
- }
- error/=dm->getSize();
- return error;
-}
-
-int main(int argc, char *argv[])
-{
- Optionpk<string> input_opt("i","input","name of the input text file");
- Optionpk<string> datum_opt("datum","datum","GPS datum of the input points","WGS84");
- Optionpk<int> s_srs_opt("s_srs","s_srs","source EPSG (integer) code of GCP input coordinates",4326);
- Optionpk<int> t_srs_opt("t_srs","t_srs","target EPSG (integer) code to project output coordinates",4326);
- Optionpk<double> focal_opt("f","focal","focal length of the camera (in m)",0.041517);
- Optionpk<double> dx_opt("dx","dx","pixel size in CCD (in micro m)",20);
- Optionpk<double> dy_opt("dy","dy","pixel size in CCD (in micro m)",20);
- Optionpk<double> x_opt("x","x","gcp x coordinate)");
- Optionpk<double> y_opt("y","y","gcp y coordinate)");
- Optionpk<double> z_opt("z","z","gcp z coordinate in m)");
- Optionpk<double> errorZ_opt("ez","ez","offset (error) in z coordinate of the GCP");
- Optionpk<int> col_opt("c","col","column of the pixel on CCD");
- Optionpk<int> row_opt("r","row","row of the pixel on CCD");
- Optionpk<double> roll_opt("roll","roll","platform attitude roll");
- Optionpk<double> pitch_opt("pitch","pitch","platform attitude pitch");
- Optionpk<double> yaw_opt("yaw","yaw","platform attitude yaw");
- Optionpk<double> xl_opt("xl","xl","platform x position in epsg:4326");
- Optionpk<double> yl_opt("yl","yl","platform y position in epsg:4326");
- Optionpk<double> zl_opt("zl","zl","platform z position (in m)");
- // Optionpk<double> bcpos_opt("bcp","bcp","bore sight offset position",0);
- Optionpk<double> bcatt_opt("bca","bca","bore sight attitude calibration offset for roll, pitch, yaw of the platform",0);
- Optionpk<double> nx_opt("nx","nx","number of columns in CCD",1441);
- Optionpk<double> ny_opt("ny","ny","number of rows in CCD",1);
- Optionpk<double> fov_opt("fov","fov","field of view (in degrees)",39.74584);
- Optionpk<double> ppx_opt("ppx","ppx","scan line principal point in x",711);
- Optionpk<double> ppy_opt("ppy","ppy","scan line principal point in y",0);
- Optionpk<char> fs_opt("fs","fs","field separator.",' ');
- Optionpk<string> output_opt("o", "output", "Output ascii file (empty: use stdout");
- Optionpk<short> sensor_opt("s", "sensor", "Sensor type (0: whiskbroom, 1: pushbroom, 2: frame, 3: apex (predefined settings))",0);
- Optionpk<double> polynome_opt("pol","polynome","coefficients for polynome to correct accross track");
- Optionpk<bool> mean_opt("m","mean","calculate mean error",false);
- Optionpk<bool> optimize_opt("opt","opt","optimize bore sight angles using GCP points in input file",false);
- Optionpk<double> lb_opt("lb","lb","lower bounds for offset bore sight angles: use -lb offset_roll -lb offset_pitch -lb offset_yaw",0);
- Optionpk<double> ub_opt("ub","ub","upper bounds for offset bore sight angles: use -ub offset_roll -ub offset_pitch -ub offset_yaw",0);
- Optionpk<unsigned int> maxit_opt("maxit","maxit","maximum number of iterations",500);
- Optionpk<double> tolerance_opt("tol","tolerance","relative tolerance for stopping criterion",0.0001);
- Optionpk<double> init_opt("is","init","initial state vector for bore sight angles: -is init_roll -is init_pitch -is init_yaw",0);
- Optionpk<double> threshold_opt("t","t","threshold for GCP error (in m)",0);
- Optionpk<bool> gcprad_opt("gcprad","gcprad","gcp coordinates",false);
- Optionpk<bool> pplrad_opt("prad","prad","platform pos coordinates",false);
- Optionpk<bool> aplrad_opt("arad","arad","platform attitude angles",false);
- Optionpk<bool> bcrad_opt("brad","brad","boresight attitude angles",false);
- Optionpk<bool> getzenith_opt("gz","getzenith","get zenith angle from platform",false);
- Optionpk<string> algorithm_opt("a", "algorithm", "optimization algorithm (see http://ab-initio.mit.edu/wiki/index.php/NLopt_Algorithms)","LN_COBYLA");
- Optionpk<short> verbose_opt("v", "verbose", "verbose mode when > 0", 0);
-
- bool doProcess;//stop process when program was invoked with help option (-h --help)
- try{
- doProcess=input_opt.retrieveOption(argc,argv);
- datum_opt.retrieveOption(argc,argv);
- s_srs_opt.retrieveOption(argc,argv);
- t_srs_opt.retrieveOption(argc,argv);
- focal_opt.retrieveOption(argc,argv);
- dx_opt.retrieveOption(argc,argv);
- dy_opt.retrieveOption(argc,argv);
- nx_opt.retrieveOption(argc,argv);
- ny_opt.retrieveOption(argc,argv);
- x_opt.retrieveOption(argc,argv);
- y_opt.retrieveOption(argc,argv);
- z_opt.retrieveOption(argc,argv);
- errorZ_opt.retrieveOption(argc,argv);
- xl_opt.retrieveOption(argc,argv);
- yl_opt.retrieveOption(argc,argv);
- zl_opt.retrieveOption(argc,argv);
- roll_opt.retrieveOption(argc,argv);
- pitch_opt.retrieveOption(argc,argv);
- yaw_opt.retrieveOption(argc,argv);
- col_opt.retrieveOption(argc,argv);
- row_opt.retrieveOption(argc,argv);
- // bcpos_opt.retrieveOption(argc,argv);
- bcatt_opt.retrieveOption(argc,argv);
- fov_opt.retrieveOption(argc,argv);
- ppx_opt.retrieveOption(argc,argv);
- ppy_opt.retrieveOption(argc,argv);
- fs_opt.retrieveOption(argc,argv);
- output_opt.retrieveOption(argc,argv);
- sensor_opt.retrieveOption(argc,argv);
- polynome_opt.retrieveOption(argc,argv);
- mean_opt.retrieveOption(argc,argv);
- optimize_opt.retrieveOption(argc,argv);
- lb_opt.retrieveOption(argc,argv);
- ub_opt.retrieveOption(argc,argv);
- maxit_opt.retrieveOption(argc,argv);
- tolerance_opt.retrieveOption(argc,argv);
- init_opt.retrieveOption(argc,argv);
- threshold_opt.retrieveOption(argc,argv);
- gcprad_opt.retrieveOption(argc,argv);
- pplrad_opt.retrieveOption(argc,argv);
- aplrad_opt.retrieveOption(argc,argv);
- bcrad_opt.retrieveOption(argc,argv);
- getzenith_opt.retrieveOption(argc,argv);
- algorithm_opt.retrieveOption(argc,argv);
- verbose_opt.retrieveOption(argc,argv);
- }
- catch(string predefinedString){
- std::cout << predefinedString << std::endl;
- exit(0);
- }
- if(!doProcess){
- std::cout << "short option -h shows basic options only, use long option --help to show all options" << std::endl;
- exit(0);//help was invoked, stop processing
- }
-
- //ID, row of GCP, column of GCP, x_GCP, y_GCP, Z_GCP, x_platform, y_platform, z_platform, platform roll, platform pitch, platform yaw of point.
- vector<int> vID;
- SensorModel::SensorModel theModel;
-
- switch(sensor_opt[0]){
- case(0):
- theModel.setModel(SensorModel::WHISKBROOM);
- break;
- case(1):
- theModel.setModel(SensorModel::PUSHBROOM);
- break;
- case(2):
- theModel.setModel(SensorModel::FRAME);
- break;
- case(3)://APEX
- theModel.setModel(SensorModel::WHISKBROOM);
- fov_opt[0]=27.984271234;
- nx_opt[0]=1000;
- ny_opt[0]=1;
- dx_opt[0]=40;
- dy_opt[0]=0;
- focal_opt[0]=0.08103;
- ppx_opt[0]=498;
- ppy_opt[0]=0;
- break;
- default:
- std::cerr << "Error: sensor type " << sensor_opt[0] << " not supported" << std::endl;
- exit(0);
- break;
- }
- theModel.setFOV(fov_opt[0]);
- theModel.setNcol(nx_opt[0]);
- theModel.setNrow(ny_opt[0]);
- theModel.setDx(dx_opt[0]);
- theModel.setDy(dy_opt[0]);
- theModel.setF(focal_opt[0]);
- theModel.setPPx(ppx_opt[0]);
- theModel.setPPy(ppy_opt[0]);
- theModel.setPolynome(polynome_opt);
- theModel.setDatum(datum_opt[0]);
- // arma::vec bc_pos(3);
- arma::vec bc_att(3);
- // while(bcpos_opt.size()<3)
- // bcpos_opt.push_back(bcpos_opt[0]);
- while(bcatt_opt.size()<3)
- bcatt_opt.push_back(bcatt_opt[0]);
- if(bcrad_opt[0]){
- // bc_pos[0]=theModel.rad2deg(bcpos_opt[0]);
- // bc_pos[1]=theModel.rad2deg(bcpos_opt[1]);
- bc_att(0)=theModel.rad2deg(bcatt_opt[0]);
- bc_att(1)=theModel.rad2deg(bcatt_opt[1]);
- bc_att(2)=theModel.rad2deg(bcatt_opt[2]);
- }
- else{
- // bc_pos(0)=bcpos_opt(0);
- // bc_pos(1)=bcpos_opt(1);
- bc_att(0)=bcatt_opt[0];
- bc_att(1)=bcatt_opt[1];
- bc_att(2)=bcatt_opt[2];
- }
- // bc_pos[2]=bcpos_opt[2];
- // theModel.setBoresightPos(bc_pos);
- theModel.setBoresightAtt(bc_att);
- if(verbose_opt[0]>1){
- // std::cout << "bore sight position offset: " << theModel.getBoresightPos() << std::endl;
- std::cout << "bore sight attitude offset: " << theModel.getBoresightAtt() << std::endl;
- }
- DataModel theDataModel(theModel);
- theDataModel.setThreshold(threshold_opt[0]);
-
- ifstream dataFile;
-
- int nrow=0;
- int ivalue=0;
- double dvalue=0;
-
- OGRSpatialReference oSourceSRS, oTargetSRS;
- OGRCoordinateTransformation *poCT;
- //input for sensor model should be universal lat lon (epsg:4326)
- oSourceSRS.importFromEPSG(s_srs_opt[0]);
- oTargetSRS.importFromEPSG(4326);
- poCT = OGRCreateCoordinateTransformation( &oSourceSRS,
- &oTargetSRS );
-
- //frame ID, row, column, GCP X, GCP Y, GCP Z, platform X, platformY, platform Z, platform roll, platform pitch, platform yaw.
- for(int ifile=0;ifile<input_opt.size();++ifile){
- if(verbose_opt[0]>1)
- std::cout << "opening file " << input_opt[ifile] << std::endl;
- dataFile.open(input_opt[ifile].c_str());
- assert(dataFile);
- if(fs_opt[0]>' '&&fs_opt[0]<='~'){//field separator is a regular character (minimum ASCII code is space, maximum ASCII code is tilde)
- string csvRecord;
- while(getline(dataFile,csvRecord)){//read a line
- istringstream csvstream(csvRecord);
- string item;
- int icol=0;//first column is id
- arma::vec thePosGCP(3);
- arma::vec thePosPlatform(3);
- arma::vec theAttPlatform(3);
- while(getline(csvstream,item,fs_opt[0])){//read a column
- if(verbose_opt[0]>1)
- std::cout << item << " ";
- switch(icol){
- case(0)://id
- ivalue=atoi(item.c_str());
- vID.push_back(ivalue);
- break;
- case(1)://row
- ivalue=atoi(item.c_str());
- theDataModel.pushRow(ivalue);
- break;
- case(2)://col
- ivalue=atoi(item.c_str());
- theDataModel.pushCol(ivalue);
- break;
- case(3)://X
- dvalue=atof(item.c_str());
- if(gcprad_opt[0])
- thePosGCP(0)=theModel.rad2deg(dvalue);
- else
- thePosGCP(0)=dvalue;
- break;
- case(4)://Y
- dvalue=atof(item.c_str());
- if(gcprad_opt[0])
- thePosGCP(1)=theModel.rad2deg(dvalue);
- else
- thePosGCP(1)=dvalue;
- break;
- case(5)://Z
- dvalue=atof(item.c_str());
- thePosGCP(2)=dvalue;
- break;
- case(6)://XL
- dvalue=atof(item.c_str());
- if(pplrad_opt[0])
- thePosPlatform(0)=theModel.rad2deg(dvalue);
- else
- thePosPlatform(0)=dvalue;
- break;
- case(7)://YL
- dvalue=atof(item.c_str());
- if(pplrad_opt[0])
- thePosPlatform(1)=theModel.rad2deg(dvalue);
- else
- thePosPlatform(1)=dvalue;
- break;
- case(8)://ZL
- dvalue=atof(item.c_str());
- thePosPlatform(2)=dvalue;
- break;
- case(9)://Roll
- dvalue=atof(item.c_str());
- if(aplrad_opt[0])
- theAttPlatform(0)=theModel.rad2deg(dvalue);
- else
- theAttPlatform(0)=dvalue;
- break;
- case(10)://Pitch
- dvalue=atof(item.c_str());
- if(aplrad_opt[0])
- theAttPlatform(1)=theModel.rad2deg(dvalue);
- else
- theAttPlatform(1)=dvalue;
- break;
- case(11)://Yaw
- dvalue=atof(item.c_str());
- if(aplrad_opt[0])
- theAttPlatform(2)=theModel.rad2deg(dvalue);
- else
- theAttPlatform(2)=dvalue;
- break;
- case(12)://track
- break;
- default:
- std::cerr << "Error: too many columns in input file " << input_opt[ifile] << std::endl;
- break;
- }
- ++icol;
- }
- if(s_srs_opt[0]!=4326){
- if(verbose_opt[0]>1)
- std::cout << "transforming: " << thePosGCP(0) << ", " << thePosGCP(1) << ", " << thePosGCP(2) << std::endl;
- poCT->Transform(1,&(thePosGCP(0)),&(thePosGCP(1)),&(thePosGCP(2)));
- if(verbose_opt[0]>1)
- std::cout << "into: " << thePosGCP(0) << ", " << thePosGCP(1) << ", " << thePosGCP(2) << std::endl;
- // poCT->Transform(1,&(thePosPlatform[0]),&(thePosPlatform[1]),&(thePosPlatform[2]));
- }
- theDataModel.pushPosGCP(thePosGCP);
- theDataModel.pushPosPlatform(thePosPlatform);
- theDataModel.pushAttPlatform(theAttPlatform);
- if(verbose_opt[0]>1)
- std::cout << endl;
- ++nrow;
- }
- }//comma separated value (or ASCII character other than space or TAB)
- else{//space or tab delimited fields
- string spaceRecord;
- while(!getline(dataFile, spaceRecord).eof()){
- if(verbose_opt[0]>1)
- std::cout << spaceRecord << std::endl;
- istringstream lineStream(spaceRecord);
- string item;
- int icol=0;
- arma::vec thePosGCP(3);
- arma::vec thePosPlatform(3);
- arma::vec theAttPlatform(3);
- while(lineStream >> item){
- if(verbose_opt[0]>1)
- std::cout << item << " ";
- istringstream itemStream(item);
- switch(icol){
- case(0)://id
- itemStream >> ivalue;
- vID.push_back(ivalue);
- break;
- case(1)://row
- itemStream >> ivalue;
- theDataModel.pushRow(ivalue);
- break;
- case(2)://col
- itemStream >> ivalue;
- theDataModel.pushCol(ivalue);
- break;
- case(3)://X
- itemStream >> dvalue;
- if(gcprad_opt[0])
- thePosGCP(0)=theModel.rad2deg(dvalue);
- else
- thePosGCP(0)=dvalue;
- break;
- case(4)://Y
- itemStream >> dvalue;
- if(gcprad_opt[0])
- thePosGCP(1)=theModel.rad2deg(dvalue);
- else
- thePosGCP(1)=dvalue;
- break;
- case(5)://Z
- itemStream >> dvalue;
- thePosGCP(2)=dvalue;
- break;
- case(6)://XL
- itemStream >> dvalue;
- if(pplrad_opt[0])
- thePosPlatform(0)=theModel.rad2deg(dvalue);
- else
- thePosPlatform(0)=dvalue;
- break;
- case(7)://YL
- itemStream >> dvalue;
- if(pplrad_opt[0])
- thePosPlatform(1)=theModel.rad2deg(dvalue);
- else
- thePosPlatform(1)=dvalue;
- break;
- case(8)://ZL
- itemStream >> dvalue;
- thePosPlatform(2)=dvalue;
- break;
- case(9)://Roll
- itemStream >> dvalue;
- if(aplrad_opt[0])
- theAttPlatform(0)=theModel.rad2deg(dvalue);
- else
- theAttPlatform(0)=dvalue;
- break;
- case(10)://Pitch
- itemStream >> dvalue;
- if(aplrad_opt[0])
- theAttPlatform(1)=theModel.rad2deg(dvalue);
- else
- theAttPlatform(1)=dvalue;
- break;
- case(11)://Yaw
- itemStream >> dvalue;
- if(aplrad_opt[0])
- theAttPlatform(2)=theModel.rad2deg(dvalue);
- else
- theAttPlatform(2)=dvalue;
- break;
- case(12)://track
- break;
- default:
- std::cerr << "Error: too many columns in input file " << input_opt[ifile] << std::endl;
- break;
- }
- ++icol;
- }
- if(s_srs_opt[0]!=4326){
- if(verbose_opt[0]>1)
- std::cout << "transforming: " << thePosGCP(0) << ", " << thePosGCP(1) << ", " << thePosGCP(2) << std::endl;
- poCT->Transform(1,&(thePosGCP(0)),&(thePosGCP(1)),&(thePosGCP(2)));
- if(verbose_opt[0]>1)
- std::cout << "into: " << thePosGCP(0) << ", " << thePosGCP(1) << ", " << thePosGCP(2) << std::endl;
- // poCT->Transform(1,&(thePosPlatform[0]),&(thePosPlatform[1]),&(thePosPlatform[2]));
- }
- theDataModel.pushPosGCP(thePosGCP);
- theDataModel.pushPosPlatform(thePosPlatform);
- theDataModel.pushAttPlatform(theAttPlatform);
- if(verbose_opt[0]>1)
- std::cout << std::endl;
- if(verbose_opt[0]>1)
- std::cout << "number of columns: " << icol << std::endl;
- ++nrow;
- }
- }
- //todo: assert sizes are all equal...
- dataFile.close();
- }//for ifile
-
- for(int ipoint=0;ipoint<xl_opt.size();++ipoint){
- if(verbose_opt[0]>1)
- std::cout << "ipoint: " << ipoint << std::endl;
- vID.push_back(ipoint);
- arma::vec thePosPlatform(3);
- arma::vec theAttPlatform(3);
- assert(row_opt.size()>ipoint);
- assert(col_opt.size()>ipoint);
- theDataModel.pushRow(row_opt[ipoint]);
- theDataModel.pushCol(col_opt[ipoint]);
- if(z_opt.size()>ipoint){
- // assert(x_opt.size()>ipoint);
- // assert(y_opt.size()>ipoint);
- arma::vec thePosGCP(3);
- if(x_opt.size()>ipoint&&y_opt.size()>ipoint){
- if(gcprad_opt[0]){
- thePosGCP(0)=theModel.rad2deg(x_opt[ipoint]);
- thePosGCP(1)=theModel.rad2deg(y_opt[ipoint]);
- }
- else{
- thePosGCP(0)=x_opt[ipoint];
- thePosGCP(1)=y_opt[ipoint];
- }
- }
- thePosGCP(2)=z_opt[ipoint];
- if(s_srs_opt[0]!=4326)
- poCT->Transform(1,&(thePosGCP(0)),&(thePosGCP(1)),&(thePosGCP(2)));
- theDataModel.pushPosGCP(thePosGCP);
- }
- assert(yl_opt.size()>ipoint);
- assert(zl_opt.size()>ipoint);
- if(pplrad_opt[0]){
- thePosPlatform(0)=theModel.rad2deg(xl_opt[ipoint]);
- thePosPlatform(1)=theModel.rad2deg(yl_opt[ipoint]);
- }
- else{
- thePosPlatform(0)=xl_opt[ipoint];
- thePosPlatform(1)=yl_opt[ipoint];
- }
- thePosPlatform(2)=zl_opt[ipoint];
- // if(s_srs_opt[0]!=4326)
- // poCT->Transform(1,&(thePosPlatform[0]),&(thePosPlatform[1]),&(thePosPlatform[2]));
- theDataModel.pushPosPlatform(thePosPlatform);
-
- assert(roll_opt.size()>ipoint);
- assert(pitch_opt.size()>ipoint);
- assert(yaw_opt.size()>ipoint);
- if(aplrad_opt[0]){
- theAttPlatform(0)=theModel.rad2deg(roll_opt[ipoint]);
- theAttPlatform(1)=theModel.rad2deg(pitch_opt[ipoint]);
- theAttPlatform(2)=theModel.rad2deg(yaw_opt[ipoint]);
- }
- else{
- theAttPlatform(0)=roll_opt[ipoint];
- theAttPlatform(1)=pitch_opt[ipoint];
- theAttPlatform(2)=yaw_opt[ipoint];
- }
- theDataModel.pushAttPlatform(theAttPlatform);
- }
-
- //todo: remove GCP points with error above threshold
- unsigned int nremoved=0;
- if(threshold_opt[0]>0){
- for(int ipoint=0;ipoint<vID.size();++ipoint){
- if(verbose_opt[0]>1)
- std::cout << "point: " << ipoint << std::endl;
- arma::vec pos_platform=theDataModel.getPosPlatform(ipoint);
- arma::vec att_platform=theDataModel.getAttPlatform(ipoint);
- arma::vec pos_calc=theDataModel.getPos(ipoint);
- arma::vec pos_gcp=theDataModel.getPosGCP(ipoint);
- ostringstream gcpss;
- double e=theDataModel.getModel().getDistGeo(pos_gcp,pos_calc);
- if(e>=theDataModel.getThreshold()){
- vID.erase(vID.begin()+ipoint);
- theDataModel.erase(ipoint);
- ++nremoved;
- }
- }
- if(!theDataModel.getSize())
- std::cerr << "Error: no data after filtering, check if input is correcly set (rad or degrees) or try to lower threshold" << std::endl;
- }
-
- if(optimize_opt[0]){
- while(lb_opt.size()<3)
- lb_opt.push_back(lb_opt[0]);
- while(ub_opt.size()<3)
- ub_opt.push_back(ub_opt[0]);
- while(init_opt.size()<3)
- init_opt.push_back(init_opt[0]);
- //todo: make nlopt::LN_COBYLA as a command line option
- //nlopt::opt opt(nlopt::LN_COBYLA,4*input_opt.size());//k,e,a,haze
- nlopt::opt optimizer=OptFactory::getOptimizer(algorithm_opt[0],3);//bore sight angles
-
- optimizer.set_min_objective(objFunction, &theDataModel);
-
- if(verbose_opt[0]>1)
- std::cout << "set lower and upper bounds" << std::endl;
- assert(lb_opt.size()==ub_opt.size());
- for(int index=0;index<lb_opt.size();++index){
- if(bcrad_opt[0]){
- lb_opt[index]=theModel.rad2deg(lb_opt[index]);
- ub_opt[index]=theModel.rad2deg(ub_opt[index]);
- }
- }
- optimizer.set_lower_bounds(lb_opt);
- optimizer.set_upper_bounds(ub_opt);
-
- if(verbose_opt[0]>1)
- std::cout << "set stopping criteria" << std::endl;
- //set stopping criteria
- if(maxit_opt[0])
- optimizer.set_maxeval(maxit_opt[0]);
- else
- optimizer.set_xtol_rel(tolerance_opt[0]);
- double minf=0;
- std::vector<double> x=init_opt;
-
- optimizer.optimize(x, minf);
- if(verbose_opt[0]){
- std::cout << "optimized with " << optimizer.get_algorithm_name() << "..." << std::endl;
- for(int index=0;index<x.size();++index){
- if(bcrad_opt[0])
- std::cout << setprecision(12) << " -bca " << theModel.deg2rad(x[index]);
- else
- std::cout << setprecision(12) << " -bca " << x[index];
- }
- std::cout << std::endl;
- }
- }
-
- double posX=0;
- double posY=0;
- ofstream outputStream;
- if(output_opt.size()){
- if(verbose_opt[0]>1)
- std::cout << "opening output file: " << output_opt[0] << std::endl;
- outputStream.open(output_opt[0].c_str());
- }
- vector<double> errorv;
-
- oSourceSRS.importFromEPSG(4326);
- oTargetSRS.importFromEPSG(t_srs_opt[0]);
- poCT = OGRCreateCoordinateTransformation( &oSourceSRS,
- &oTargetSRS );
- for(int ipoint=0;ipoint<vID.size();++ipoint){
- if(verbose_opt[0]>1)
- std::cout << "point: " << ipoint << std::endl;
- arma::vec pos_platform=theDataModel.getPosPlatform(ipoint);
- arma::vec att_platform=theDataModel.getAttPlatform(ipoint);
-
- arma::vec pos_calc;
- arma::vec pos_gcp;
- if(input_opt.size()||x_opt.size())
- pos_gcp=theDataModel.getPosGCP(ipoint);
- else
- pos_gcp=theDataModel.getPos(ipoint);
- if(errorZ_opt.size())
- pos_calc=theDataModel.getModel().getPos(pos_platform,att_platform,theDataModel.getRow(ipoint),theDataModel.getCol(ipoint),theDataModel.getHeight(ipoint)+errorZ_opt[0]);
- else
- pos_calc=theDataModel.getPos(ipoint);
- double e=theDataModel.getModel().getDistGeo(pos_gcp,pos_calc);
- if(theDataModel.getThreshold()&&e>theDataModel.getThreshold())
- continue;
- errorv.push_back(e);
-
- ostringstream gcpss;
- gcpss.precision(12);
- poCT->Transform(1,&(pos_gcp(0)),&(pos_gcp(1)),&(pos_gcp(2)));
- gcpss << pos_gcp(0) << " " << pos_gcp(1) << " " << pos_gcp(2) << " ";
- gcpss << errorv.back() << " ";
-
- poCT->Transform(1,&(pos_calc(0)),&(pos_calc(1)),&(pos_calc(2)));
-
- if(output_opt.size()){
- outputStream << setprecision(12) << vID[ipoint] << " " << theDataModel.getRow(ipoint) << " " << theDataModel.getCol(ipoint) << " " << pos_calc(0) << " " << pos_calc(1) << " " << pos_calc(2) << " " << gcpss.str();
- if(getzenith_opt[0])
- outputStream << " " << theDataModel.getModel().getZenith(att_platform,theDataModel.getRow(ipoint),theDataModel.getCol(ipoint));
- outputStream << std::endl;
- }
- else{
- std::cout << setprecision(12) << vID[ipoint] << " " << theDataModel.getRow(ipoint) << " " << theDataModel.getCol(ipoint) << " " << pos_calc(0) << " " << pos_calc(1) << " " << pos_calc(2) << " " << gcpss.str();
- if(getzenith_opt[0])
- std::cout << " " << theDataModel.getModel().getZenith(att_platform,theDataModel.getRow(ipoint),theDataModel.getCol(ipoint));
- std::cout << std::endl;
- }
- }
- if(verbose_opt[0]){
- if(output_opt.size())
- outputStream << "Number of GCP above threshold removed: " << nremoved << std::endl;
- else
- std::cout << "Number of GCP above threshold removed: " << nremoved << std::endl;
- }
- if(mean_opt[0]){
- double theMean=0;
- double theVar=0;
- statfactory::StatFactory stat;
- stat.meanVar(errorv,theMean,theVar);
- if(output_opt.size())
- outputStream << setprecision(12) << "mean stddev: " << theMean << " " << sqrt(theVar) << std::endl;
- else
- std::cout << setprecision(12) << "mean stdev: " << theMean << " " << sqrt(theVar) << std::endl;
- }
- if(output_opt.size())
- outputStream.close();
-}
diff --git a/src/apps/pksetchandelier.cc b/src/apps/pksetchandelier.cc
deleted file mode 100644
index 95bca22..0000000
--- a/src/apps/pksetchandelier.cc
+++ /dev/null
@@ -1,224 +0,0 @@
-/**********************************************************************
-pksetchandelier.cc: program to apply model parameters for brdf correction
-Copyright (C) 2008-2012 Pieter Kempeneers
-
-This file is part of pktools
-
-pktools is free software: you can redistribute it and/or modify
-it under the terms of the GNU General Public License as published by
-the Free Software Foundation, either version 3 of the License, or
-(at your option) any later version.
-
-pktools is distributed in the hope that it will be useful,
-but WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-GNU General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with pktools. If not, see <http://www.gnu.org/licenses/>.
-***********************************************************************/
-#include <assert.h>
-#include <string>
-#include <list>
-#include <iostream>
-#include "base/PointData.h"
-#include "imageclasses/ImgReaderOgr.h"
-#include "imageclasses/ImgWriterOgr.h"
-#include "Optionpk.h"
-
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
-
-int main(int argc, char *argv[])
-{
- std::string versionString="version ";
- versionString+=VERSION;
- versionString+=", Copyright (C) 2008-2012 Pieter Kempeneers.\n\
- This program comes with ABSOLUTELY NO WARRANTY; for details type use option -h.\n\
- This is free software, and you are welcome to redistribute it\n\
- under certain conditions; use option --license for details.";
- Optionpk<bool> version_opt("\0","version",versionString,false);
- Optionpk<bool> license_opt("lic","license","show license information",false);
- Optionpk<bool> help_opt("h","help","shows this help info",false);
- Optionpk<bool> todo_opt("\0","todo","",false);
- Optionpk<string> input_opt("i", "input", "Reflectance input","");
- Optionpk<string> geom_opt("g", "geom", "Geometry (ogr vector) file","");
- Optionpk<string> output_opt("o", "output", "corrected image output","");
- Optionpk<string> otype_opt("ot", "otype", "Data type for output image ({Byte/Int16/UInt16/UInt32/Int32/Float32/Float64/CInt16/CInt32/CFloat32/CFloat64}). Empty string: inherit type from input image", "Float32");
- Optionpk<string> oformat_opt("of", "oformat", "Output image format (see also gdal_translate). Empty string: inherit from input image", "GTiff");
- Optionpk<string> option_opt("co", "co", "options: NAME=VALUE [-co COMPRESS=LZW] [-co INTERLEAVE=BAND]");
- Optionpk<unsigned short> band_opt("b", "band", "field name of Reflectance band in reflectance file",0);
- Optionpk<unsigned short> sza_opt("sza", "sza", "band number (starting from 0) for Sun Zenith Angle in geometry image file",3);
- // Optionpk<unsigned short> vza_opt("vza", "vza", "band number (starting from 0) for Sun Zenith Angle in geometry image file",5);
- Optionpk<unsigned short> saa_opt("saa", "saa", "band number (starting from 0) for Sun Azimuth Angle in geometry image file",6);
- Optionpk<unsigned short> vaa_opt("vaa", "vaa", "band number (starting from 0) for View Azimuth Angle in geometry image file",7);
- Optionpk<double> tau_opt("tau", "tau", "Optical depth",0.2);
- Optionpk<double> deltaT_opt("deltaT", "deltaT", "Exposure time",1.0);
- Optionpk<double> state_opt("s","state","state vector: -s k -s e -s a -s haze",0);
- Optionpk<bool> model_opt("m","model","forward model reflectance",false);
- Optionpk<double> flag_opt("f", "flag", "No value flag",1000);
- Optionpk<short> verbose_opt("v", "verbose", "verbose", 0);
-
- version_opt.retrieveOption(argc,argv);
- license_opt.retrieveOption(argc,argv);
- help_opt.retrieveOption(argc,argv);
- todo_opt.retrieveOption(argc,argv);
-
- input_opt.retrieveOption(argc,argv);
- geom_opt.retrieveOption(argc,argv);
- output_opt.retrieveOption(argc,argv);
- otype_opt.retrieveOption(argc,argv);
- oformat_opt.retrieveOption(argc,argv);
- option_opt.retrieveOption(argc,argv);
- band_opt.retrieveOption(argc,argv);
- sza_opt.retrieveOption(argc,argv);
- // vza_opt.retrieveOption(argc,argv);
- saa_opt.retrieveOption(argc,argv);
- vaa_opt.retrieveOption(argc,argv);
- tau_opt.retrieveOption(argc,argv);
- deltaT_opt.retrieveOption(argc,argv);
- state_opt.retrieveOption(argc,argv);
- model_opt.retrieveOption(argc,argv);
- flag_opt.retrieveOption(argc,argv);
- verbose_opt.retrieveOption(argc,argv);
-
- if(version_opt[0]||todo_opt[0]){
- cout << version_opt.getHelp() << endl;
- cout << "todo: " << todo_opt.getHelp() << endl;
- exit(0);
- }
- if(license_opt[0]){
- cout << Optionpk<bool>::getGPLv3License() << endl;
- exit(0);
- }
- if(help_opt[0]){
- cout << "usage: pkbrdf -i input -o output -g geom [OPTIONS]" << endl;
- exit(0);
- }
-
- if(verbose_opt[0])
- std::cout << "state_opt.size(): " << state_opt.size() << std::endl;
- assert(state_opt.size()==4);
- ImgReaderGdal inputReader(input_opt[0]);
- ImgReaderGdal geomReader(geom_opt[0]);
- vector<double> inputBuffer(inputReader.nrOfCol());
- vector<double> szaBuffer(geomReader.nrOfCol());
- vector<double> saaBuffer(geomReader.nrOfCol());
- // vector<double> vzaBuffer(geomReader.nrOfCol());
- vector<double> vaaBuffer(geomReader.nrOfCol());
- vector<double> fiBuffer(geomReader.nrOfCol());
- vector<double> xBuffer(inputReader.nrOfCol());//x pos
- vector<double> yBuffer(inputReader.nrOfCol());//y pos
- ImgWriterGdal outputWriter;
- string imageType=inputReader.getImageType();
- if(oformat_opt[0]!="")//default
- imageType=oformat_opt[0];
- GDALDataType theType=GDT_Unknown;
- if(verbose_opt[0]){
- std::cout << "Image type: " << imageType << std::endl;
- std::cout << "possible output data types: ";
- }
- for(int iType = 0; iType < GDT_TypeCount; ++iType){
- if(verbose_opt[0])
- cout << " " << GDALGetDataTypeName((GDALDataType)iType);
- if( GDALGetDataTypeName((GDALDataType)iType) != NULL
- && EQUAL(GDALGetDataTypeName((GDALDataType)iType),
- otype_opt[0].c_str()))
- theType=(GDALDataType) iType;
- }
- if(theType==GDT_Unknown)
- theType=inputReader.getDataType();
-
- if(verbose_opt[0]){
- std::cout << std::endl << "Output data type: " << GDALGetDataTypeName(theType) << std::endl;
- std::cout << "opening output image for writing: " << output_opt[0] << std::endl;
- }
- if(option_opt.findSubstring("INTERLEAVE=")==option_opt.end()){
- string theInterleave="INTERLEAVE=";
- theInterleave+=inputReader.getInterleave();
- option_opt.push_back(theInterleave);
- }
- try{
- outputWriter.open(output_opt[0],inputReader.nrOfCol(),inputReader.nrOfRow(),1,theType,imageType,option_opt);
- outputWriter.setProjection(inputReader.getProjection());
- outputWriter.copyGeoTransform(inputReader);
- }
- catch(string errorstring){
- cout << errorstring << endl;
- exit(1);
- }
- vector<double> writeBuffer(outputWriter.nrOfCol());
-
- PointData::m_tau=tau_opt[0];
- PointData::m_deltaT=deltaT_opt[0];
-
- const char* pszMessage;
- void* pProgressArg=NULL;
- GDALProgressFunc pfnProgress=GDALTermProgress;
- double progress=0;
- pfnProgress(progress,pszMessage,pProgressArg);
- for(int irow=0;irow<inputReader.nrOfRow();++irow){
- try{
- geomReader.readData(szaBuffer,GDT_Float64,irow,sza_opt[0]);
- geomReader.readData(saaBuffer,GDT_Float64,irow,saa_opt[0]);
- // geomReader.readData(vzaBuffer,GDT_Float64,irow,vza_opt[0]);
- geomReader.readData(vaaBuffer,GDT_Float64,irow,vaa_opt[0]);
- }
- catch(string errorstring){
- cout << errorstring << endl;
- exit(1);
- }
- //convert angles to radiance values
- for(int icol=0;icol<inputReader.nrOfCol();++icol){
- inputReader.image2geo(icol,irow,xBuffer[icol],yBuffer[icol]);
- double theRelativeAzimuth=saaBuffer[icol]-vaaBuffer[icol];
- fiBuffer[icol]=theRelativeAzimuth;
- }
- assert(band_opt[0]<inputReader.nrOfBand());
- try{
- inputReader.readData(inputBuffer,GDT_Float64,irow,band_opt[0]);
- }
- catch(string errorstring){
- cout << errorstring << endl;
- exit(1);
- }
- PointData pd;
- State theState;
- for(int index=0;index<inputBuffer.size();++index){
- pd.setReflectance(inputBuffer[index]);
- pd.setCosSZA(szaBuffer[index]);
- pd.setCosFi(fiBuffer[index]);
- double centreX,centreY;
- inputReader.getCentrePos(centreX,centreY);
- double ulx=inputReader.getUlx();
- double uly=inputReader.getUly();
- double r0=((ulx-centreX)*(ulx-centreX)+(uly-centreY)*(uly-centreY));
- double r=((xBuffer[index]-centreX)*(xBuffer[index]-centreX)+(yBuffer[index]-centreY)*(yBuffer[index]-centreY));
- pd.setR(r/r0);
- theState.k=state_opt[0];
- theState.e=state_opt[1];
- theState.a=state_opt[2];
- theState.haze=state_opt[3];
- if(model_opt[0])
- writeBuffer[index]=pd.modelReflectance(theState);
- else
- writeBuffer[index]=pd.correctReflectance(theState);
- }
- try{
- outputWriter.writeData(writeBuffer,GDT_Float64,irow,0);
- }
- catch(string errorstring){
- cout << errorstring << endl;
- exit(1);
- }
- progress=(1.0+irow)/inputReader.nrOfRow();
- assert(progress>=0);
- assert(progress<=1);
- pfnProgress(progress,pszMessage,pProgressArg);
- }
- inputReader.close();
- geomReader.close();
- if(output_opt[0]!="")
- outputWriter.close();
-}
diff --git a/src/base/Makefile.am b/src/base/Makefile.am
index 927e830..25c6dfa 100644
--- a/src/base/Makefile.am
+++ b/src/base/Makefile.am
@@ -20,10 +20,9 @@ lib_LTLIBRARIES = libbase.la
libbase_ladir = $(includedir)/base
# the list of header files that belong to the library (to be installed later)
-libbase_la_HEADERS = $(top_srcdir)/config.h PointData.h Vector2d.h IndexValue.h Optionpk.h PosValue.h
+libbase_la_HEADERS = $(top_srcdir)/config.h Vector2d.h IndexValue.h Optionpk.h PosValue.h
# the sources to add to the library and to add to the source distribution
-libbase_la_SOURCES = $(libbase_la_HEADERS) PointData.cc
###############################################################################
# list of sources for the binaries
@@ -33,5 +32,5 @@ pktestOption_SOURCES = pktestOption.cc
# HEADERS USED FOR BUILDING BUT NOT TO BE INSTALLED
###############################################################################
-#noinst_HEADERS = Optionpk.h PosValue.h Vector2d.h PointData.h
+#noinst_HEADERS = Optionpk.h PosValue.h Vector2d.h
###############################################################################
diff --git a/src/base/PointData.cc b/src/base/PointData.cc
deleted file mode 100644
index 3d55767..0000000
--- a/src/base/PointData.cc
+++ /dev/null
@@ -1,34 +0,0 @@
-#include "PointData.h"
-
-double PointData::m_tau=0;
-double PointData::m_deltaT=0;
-double PointData::m_epsilon=0.00001;
-double PointData::m_flag=1000;
-double PointData::m_residual=0;
-
-PointData::PointData(){
-}
-
-PointData::~PointData(){
-}
-
-double PointData::correctReflectance(const State &state){
- double ecosfi=state.e*m_cosfi;
- double optical_depth=(m_cosza*m_cosza>m_epsilon*m_epsilon) ? exp(-m_tau/m_cosza) : 0;
- double hotspot=(ecosfi<1) ? (state.k+(1-state.k)*(1-state.e*m_cosza)/(1-ecosfi)) : 1;
- double radial=(1+state.a*m_r);
- double correctedReflectance;
- correctedReflectance=(m_reflectance-state.haze)/m_cosza/optical_depth/m_deltaT/hotspot/radial;
- return correctedReflectance;
- // return (m_reflectance-state.haze)/((m_cosza*exp(-m_tau/m_cosza)*m_deltaT*(state.k+(1-state.k)*(1-state.e*m_cosza)/(1-ecosfi))*(1+state.a*m_r)));
-}
-
-double PointData::modelReflectance(const State &state){
- double ecosfi=state.e*m_cosfi;
- double optical_depth=exp(-m_tau/m_cosza);
- double hotspot=state.k+(1-state.k)*(1-state.e*m_cosza)/(1-ecosfi);
- double radial=(1+state.a*m_r);
- return m_reflectance*m_cosza*optical_depth*m_deltaT*hotspot*radial+state.haze;
- // return (m_reflectance-state.haze)/((m_cosza*exp(-m_tau/m_cosza)*m_deltaT*(state.k+(1-state.k)*(1-state.e*m_cosza)/(1-ecosfi))*(1+state.a*m_r)));
-}
-
diff --git a/src/base/PointData.h b/src/base/PointData.h
deleted file mode 100644
index 08c5a1a..0000000
--- a/src/base/PointData.h
+++ /dev/null
@@ -1,51 +0,0 @@
-#ifndef _POINTDATA_
-#define _POINTDATA_
-#include <assert.h>
-#include <math.h>
-#include <vector>
-
-struct State{
- double k;
- double e;
- double a;
- double haze;
-};
-
-class PointData{
- public:
- PointData();//constructor
- ~PointData();//destructor
- void setCosFi(double fi){
- if(fi>180)
- fi-=360;
- else if(fi<-180)
- fi+=360;
- m_cosfi=cos(deg2rad(fi));
- };
- void setCosSZA(double sza){
- assert(sza>=0);
- m_cosza=cos(deg2rad(sza));
- };
- double getReflectance() const{return m_reflectance;};
- void setReflectance(double reflectance) {m_reflectance=reflectance;};
- double getImage() const{return m_image;};
- void setImage(int image) {m_image=image;};
- void setR(double r){m_r=r;};
- double getR() const{return m_r;};
- double correctReflectance(const State &state);
- double modelReflectance(const State &state);
- static double m_deltaT;
- static double m_tau;
- static double m_epsilon;
- static double m_flag;
- static double m_residual;
- private:
- double deg2rad(double angle) const{return acos(-1)*angle/180.0;};
- double m_reflectance;
- double m_cosza;
- double m_cosfi;
- double m_r;
- int m_image;
-};
-
-#endif //_POINTDATA_
diff --git a/src/models/Makefile.am b/src/models/Makefile.am
deleted file mode 100644
index 36cb684..0000000
--- a/src/models/Makefile.am
+++ /dev/null
@@ -1,40 +0,0 @@
-AM_LDFLAGS = $(GSL_LIBS) $(GDAL_LDFLAGS) @AM_LDFLAGS@
-AM_CXXFLAGS = -I$(top_srcdir)/src $(GDAL_CFLAGS) @AM_CXXFLAGS@
-
-###############################################################################
-# THE PROGRAMS TO BUILD
-###############################################################################
-
-# the program to build (the names of the final binaries)
-#do not want to install pktestoption
-noinst_PROGRAMS = pktestProspect
-
-###############################################################################
-# THE LIBRARIES TO BUILD
-###############################################################################
-
-# the library names to build (note we are building static libs only)
-#noinst_LIBRARIES = libmodels.a
-lib_LTLIBRARIES = libprospect.la libmodels.la
-
-# where to install the headers on the system
-libmodels_ladir = $(includedir)/models
-
-# the list of header files that belong to the library (to be installed later)
-libmodels_la_HEADERS = Prospect.h
-
-if USE_NLOPT
-libmodels_la_HEADERS += SensorModel.h
-endif
-
-# the sources to add to the library and to add to the source distribution
-libmodels_la_SOURCES = $(libmodels_la_HEADERS) Prospect.cc
-libmodels_la_LIBADD = $(top_builddir)/src/models/libprospect.la
-
-libprospect_la_SOURCES = dataSpec_P5B.f90 prospect_5B.f90 tav_abs.f90
-###############################################################################
-
-# list of sources for the binaries
-pktestProspect_SOURCES = pktestProspect.cc
-#pktestProspect_LDADD = $(top_builddir)/src/models/libmodels.la $(top_builddir)/src/models/libprospect.a libf95.a
-pktestProspect_LDADD = $(top_builddir)/src/models/libmodels.la $(top_builddir)/src/models/libprospect.la
diff --git a/src/models/Prospect.cc b/src/models/Prospect.cc
deleted file mode 100644
index 155359f..0000000
--- a/src/models/Prospect.cc
+++ /dev/null
@@ -1,31 +0,0 @@
-/**********************************************************************
-Prospect.cc: class for radiative transfer leaf PROSPECT
-Copyright (C) 2008-2013 Pieter Kempeneers
-
-This file is part of pktools
-
-pktools is free software: you can redistribute it and/or modify
-it under the terms of the GNU General Public License as published by
-the Free Software Foundation, either version 3 of the License, or
-(at your option) any later version.
-
-pktools is distributed in the hope that it will be useful,
-but WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-GNU General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with pktools. If not, see <http://www.gnu.org/licenses/>.
-***********************************************************************/
-#include "Prospect.h"
-
-void prospect::Prospect::getLeafSpectrum(std::vector<double>& leafReflectance, std::vector<double>& leafTransmittance) const{
- double RT[2][2101];//400 nm - 2500 nm
- prospect_5b_(m_N,m_Cab,m_Car,m_Cbrown,m_Cw,m_Cm,RT);
- leafReflectance.resize(2101);
- leafTransmittance.resize(2101);
- for(int index=0;index<2101;++index){
- leafReflectance[index]=RT[0][index];
- leafTransmittance[index]=RT[1][index];
- }
-}
diff --git a/src/models/Prospect.h b/src/models/Prospect.h
deleted file mode 100644
index f92f5ff..0000000
--- a/src/models/Prospect.h
+++ /dev/null
@@ -1,58 +0,0 @@
-/**********************************************************************
-Prospect.h: class for radiative transfer leaf PROSPECT
-Copyright (C) 2008-2013 Pieter Kempeneers
-
-This file is part of pktools
-
-pktools is free software: you can redistribute it and/or modify
-it under the terms of the GNU General Public License as published by
-the Free Software Foundation, either version 3 of the License, or
-(at your option) any later version.
-
-pktools is distributed in the hope that it will be useful,
-but WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-GNU General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with pktools. If not, see <http://www.gnu.org/licenses/>.
-***********************************************************************/
-#ifndef _PROSPECT_
-#define _PROSPECT_
-#include <assert.h>
-#include <math.h>
-#include <vector>
-
-extern "C" {
- void prospect_5b_(const double& N,const double& Cab,const double& Car, const double& Cbrown, const double& Cw, const double& Cm, double RT[2][2101]);
-}
-
-namespace prospect
-{
- class Prospect{
- public:
- Prospect(void){};
- ~Prospect(void){};
- void setN(double n){m_N=n;};
- void setCab(double cab){m_Cab=cab;};
- void setCar(double car){m_Car=car;};
- void setCbrown(double cbrown){m_Cbrown=cbrown;};
- void setCw(double cw){m_Cw=cw;};
- void setCm(double cm){m_Cm=cm;};
- double getN() const {return m_N;};
- double getCab() const {return m_Cab;};
- double getCar() const {return m_Car;};
- double getCbrown() const {return m_Cbrown;};
- double getCw() const {return m_Cw;};
- double getCm() const {return m_Cm;};
- void getLeafSpectrum(std::vector<double>& leafReflectance,std::vector<double>& leafTransmittance) const;
- private:
- double m_N; // structure coefficient
- double m_Cab; // chlorophyll content (g.cm-2)
- double m_Car; // carotenoid content (g.cm-2)
- double m_Cbrown; // brown pigment content (arbitrary units)
- double m_Cw; // EWT (cm)
- double m_Cm; // LMA (g.cm-2)
- };
-}
-#endif //_PROSPECT_
diff --git a/src/models/SensorModel.h b/src/models/SensorModel.h
deleted file mode 100644
index 77209e2..0000000
--- a/src/models/SensorModel.h
+++ /dev/null
@@ -1,399 +0,0 @@
-/**********************************************************************
-SensorModel.h: class for sensor model (calculate position on earth from pos and attidude information)
-Copyright (C) 2008-2012 Pieter Kempeneers
-
-This file is part of pktools
-
-pktools is free software: you can redistribute it and/or modify
-it under the terms of the GNU General Public License as published by
-the Free Software Foundation, either version 3 of the License, or
-(at your option) any later version.
-
-pktools is distributed in the hope that it will be useful,
-but WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-GNU General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with pktools. If not, see <http://www.gnu.org/licenses/>.
-***********************************************************************/
-#ifndef _SENSORMODEL_
-#define _SENSORMODEL_
-#include <assert.h>
-#include <math.h>
-#include <vector>
-#include <iostream>
-#include <armadillo>
-#include "ogr_spatialref.h"
-
-#ifndef PI
-#define PI 3.1415926535897932384626433832795
-#endif
-
-// #ifndef DEG2RAD
-// #define DEG2RAD(DEG) (DEG/180.0*PI)
-// #endif
-
-namespace SensorModel
-{
- enum Type { PUSHBROOM=0, WHISKBROOM=1 , FRAME=2};
-
- class SensorModel{
- public:
- SensorModel(void){
- m_bcpos.set_size(3);
- m_bcatt.set_size(3);
- m_bcpos[0]=0;
- m_bcpos[1]=0;
- m_bcpos[2]=0;
- m_bcatt[0]=0;
- m_bcatt[1]=0;
- m_bcatt[2]=0;
- };
- SensorModel(int theModel) : m_model(theModel){
- m_bcpos.set_size(3);
- m_bcatt.set_size(3);
- m_bcpos[0]=0;
- m_bcpos[1]=0;
- m_bcpos[2]=0;
- m_bcatt[0]=0;
- m_bcatt[1]=0;
- m_bcatt[2]=0;
- };
- ~SensorModel(){};
- // double deg2rad(double angle) const{return acos(-1)*angle/180.0;};
- static double deg2rad(double angle) {return PI*angle/180.0;};
- static double rad2deg(double angle) {return 180*angle/PI;};
- void setModel(int theModel){m_model=theModel;};
- int getModel() const {return m_model;};
- void setFOV(double fovDEG){m_fov=deg2rad(fovDEG);};
- void setNcol(int ncol){m_ncol=ncol;};
- void setNrow(int nrow){m_nrow=nrow;};
- void setDx(double dx){m_dx=dx;};
- void setDy(double dy){m_dy=dy;};
- void setF(double fc){m_fc=fc;};
- double getF() const {return m_fc;};
- void setPPx(double ppx){m_ppx=ppx;};
- void setPPy(double ppy){m_ppy=ppy;};
- void setBoresightPos(const arma::vec& bcpos){m_bcpos=bcpos;};
- void setBoresightAtt(const arma::vec& bcatt){m_bcatt=bcatt;};
- arma::vec getBoresightPos() const{return m_bcpos;};
- arma::vec getBoresightAtt() const{return m_bcatt;};
- void setPolynome(const std::vector<double>& polynome) {m_polynome=polynome;};
- void setDatum(const std::string& theDatum="WGS84"){m_spatialref.SetWellKnownGeogCS(theDatum.c_str());};
- double getZenith(const arma::vec& att_platform, int row, int column) const{
- arma::vec normallevel(3);
- arma::vec normalplatform(3);
-
- arma::vec apl_deg=att_platform;
- apl_deg+=m_bcatt;
- arma::vec apl_rad(3);
- apl_rad(0)=deg2rad(apl_deg(0));//roll
- apl_rad(1)=deg2rad(apl_deg(1));//pitch
- apl_rad(2)=deg2rad(apl_deg(2));//yaw
-
- if(getModel()==PUSHBROOM){
- arma::vec scanAngle(2);
- scanAngle=scanAngle_PB(row,column);
- apl_rad(0)+=scanAngle(1);
- apl_rad(1)+=scanAngle(0);
- }
- else if(getModel()==WHISKBROOM){
- apl_rad(0)+=0;
- apl_rad(1)=tan(scanAngle_WB(column));
- }
- else//not implemented?
- assert(0);
- normallevel(0)=0;
- normallevel(1)=0;
- normallevel(2)=-1;
- arma::Mat<double> rotM(3,3);
- rotM=getRz(apl_rad(2))*getRy(apl_rad(1))*getRx(apl_rad(0));
- normalplatform=rotM*normallevel;
- return rad2deg(acos(-normalplatform(2)));
- };
- arma::vec getPos(const arma::vec& pos_platform, const arma::vec& att_platform, int row, int column, double elevation) const{
- arma::vec thePosition(3);
- arma::vec pos_ellips(3);
- arma::vec ppl_deg=pos_platform;
- arma::vec apl_deg=att_platform;
- ppl_deg+=m_bcpos;
- apl_deg+=m_bcatt;
- arma::vec ppl_rad(3);
- arma::vec apl_rad(3);
- ppl_rad(0)=deg2rad(ppl_deg(0));
- ppl_rad(1)=deg2rad(ppl_deg(1));
- ppl_rad(2)=ppl_deg(2);//add geoid elevation if necessary...
- apl_rad(0)=deg2rad(apl_deg(0));//roll
- apl_rad(1)=deg2rad(apl_deg(1));//pitch
- apl_rad(2)=deg2rad(apl_deg(2));//yaw
- arma::vec pos_ecef=pECEF(ppl_rad,apl_rad,row,column);
- pos_ellips=ecef2geo(pos_ecef);
- pos_ellips(2)=0;
- thePosition=getXYatZ(elevation,ppl_deg,pos_ellips);
- return thePosition;
- };
- double getDistGeo(const arma::vec& pos1, const arma::vec& pos2) const{
- double lon1=pos1(0);
- double lat1=pos1(1);
- double lon2=pos2(0);
- double lat2=pos2(1);
- double result;
- //simplified formula (spherical approximation)
- // result=2*asin(sqrt(pow(sin((lat1-lat2)/2),2) + cos(lat1)*cos(lat2)*pow(sin((lon1-lon2)/2),2)));
- //using loxodromes
- result=(pos1(1)==pos2(1)) ? lox2(deg2rad(pos1(1)),deg2rad(pos1(0)),deg2rad(pos2(0))) : lox1(deg2rad(pos1(1)),deg2rad(pos1(0)),deg2rad(pos2(1)),deg2rad(pos2(0)));
- return result;
- };
- arma::vec geo2ecef(const arma::vec& pos_geo) const{
- arma::vec pos_ecef(3);
- double nu=getNu(pos_geo(1));
- double f1=(nu+pos_geo(2))*cos(pos_geo(1));
- double e1=getE1();
- pos_ecef(0)=f1*cos(pos_geo(0));
- pos_ecef(1)=f1*sin(pos_geo(0));
- pos_ecef(2)=(nu*(1-e1*e1)+pos_geo(2))*sin(pos_geo(1));
- return pos_ecef;
- };
- arma::vec ecef2geo(const arma::vec& pos_ecef) const{
- arma::vec pos_geo(3);
- pos_geo(0)=atan2(pos_ecef(1),pos_ecef(0));
- while(pos_geo(0)>2*PI)
- pos_geo(0)-=2*PI;
- double r=sqrt(pos_ecef(0)*pos_ecef(0)+pos_ecef(1)*pos_ecef(1));
- double f_earth=1.0/m_spatialref.GetInvFlattening();
- double e2=f_earth*(2-f_earth);
- double Ae=m_spatialref.GetSemiMajor();
- pos_geo(1)=atan(pos_ecef(2)/r);
- double c=1;
- double lat=1E+30;
- int iterations=0;
- while(sqrt((pos_geo(1)-lat)*(pos_geo(1)-lat))>1E-15){
- ++iterations;
- lat=pos_geo(1);
- c=1.0/sqrt(1-e2*sin(lat)*sin(lat));
- pos_geo(1)=atan((pos_ecef(2)+Ae*c*e2*sin(lat))/r);
- }
- pos_geo(2)=r/cos(pos_geo(1))-Ae*c;
- pos_geo(0)=rad2deg(pos_geo(0));
- pos_geo(1)=rad2deg(pos_geo(1));
- return pos_geo;
- };
- private:
- //line function; In this function, point 0 is the platform position, point 1 is the ellipsoid position
- arma::vec getXYatZ (double inZ, arma::vec p0, arma::vec p1) const{
- arma::vec posatz(3);
- posatz(0)=p0(0)+(p0(0)-p1(0))/(p0(2)-p1(2))*(inZ-p0(2));
- posatz(1)=p0(1)+(p0(1)-p1(1))/(p0(2)-p1(2))*(inZ-p0(2));
- posatz(2)=inZ;
- return posatz;
- };
- //get First eccentricity of the Earth's surface
- double getE1() const{double f_earth=1.0/m_spatialref.GetInvFlattening();return sqrt(f_earth*(2-f_earth));};
- //get Second eccentricity of the Earth's surface
- double getEearth() const{double Ae=m_spatialref.GetSemiMajor(); double Be=m_spatialref.GetSemiMinor();return sqrt(Ae*Ae-Be*Be)/Ae;};
- //get Radius of curvature normal to the meridian
- double getNu(double lat) const{double sinlat=sin(lat);return m_spatialref.GetSemiMajor()/sqrt(1-getE1()*getE1()*sinlat*sinlat);};
- //get Ellipsoid radius in the plane of the meridian
- double getRo(double lat) const{double e1=getE1();double sinlat=sin(lat);return (m_spatialref.GetSemiMajor()*(1-e1*e1))/sqrt(1-getE1()*getE1()*sinlat*sinlat);};
- double getEccentricity() const{return sqrt(1-m_spatialref.GetSemiMinor()*m_spatialref.GetSemiMinor()/m_spatialref.GetSemiMajor()/m_spatialref.GetSemiMajor());};
- //Loxodromes are closely related with Mercator projection as the main feature of this projection is that loxodromes are projected as straight lines. To find the angle of the loxodrome between two points one has simply to project the two points using a Mercator projection and calculate the angle of the line joining the two resulting points. Formulas for Mercator projection on the ellipsoid can be found in many places (for instance "Map Projections - A Working Manual" by J.P. Snyder, U [...]
- //project lat lon to Mercator projection
- void getMxy(double& lon, double& lat) const{
- // lon=m_spatialref.GetSemiMajor()*lon;
- // double sinlat=sin(lat);
- // double f1=(1+sinlat)/(1-sinlat);
- // double eps=getEccentricity();
- // double f2=pow((1-eps*sinlat)/(1+eps*sinlat),eps);
- // lat=m_spatialref.GetSemiMajor()/2.0*log(f1*f2);
- // //test
- // std::cout << "Jan Mercator lat, lon: " << lat << ", " << lon << std::endl;
- OGRSpatialReference oSourceSRS, oTargetSRS;
- OGRCoordinateTransformation *poCT;
- oSourceSRS.importFromEPSG(4326);
- oTargetSRS.importFromEPSG(3395);
- poCT = OGRCreateCoordinateTransformation( &oSourceSRS,
- &oTargetSRS );
- lon=rad2deg(lon);
- lat=rad2deg(lat);
- poCT->Transform(1,&lon,&lat);
- };
-
- arma::vec pECEF(const arma::vec& pos, const arma::vec& attitude, int row, int column) const{
- arma::vec A(3);
- arma::Mat<double> B(3,3);
- B(0,0)=-sin(pos(1))*cos(pos(0));
- B(0,1)=-sin(pos(0));
- B(0,2)=-cos(pos(1))*cos(pos(0));
- B(1,0)=-sin(pos(1))*sin(pos(0));
- B(1,1)=cos(pos(0));
- B(1,2)=-cos(pos(1))*sin(pos(0));
- B(2,0)=cos(pos(1));
- B(2,1)=0;
- B(2,2)=-sin(pos(1));
- A=geo2ecef(pos);
- arma::vec C(3);
- A+=B*SV(row,column,attitude,pos(2));
- return A;
- };
- arma::vec getIV(int row, int column) const{
- arma::vec iv(3);
- iv(0)=0;
- arma::vec scanAngle(2);
- if(getModel()==PUSHBROOM){
- scanAngle=scanAngle_PB(row,column);
- iv(0)=scanAngle(1);
- iv(1)=tan(scanAngle(0));
- }
- else if(getModel()==WHISKBROOM){
- iv(0)=0;
- iv(1)=tan(scanAngle_WB(column));
- }
- else//not implemented?
- assert(0);
- iv(2)=-1.0;
- return iv;
- };
- double scanAngle_WB(int column) const{
- return (-m_fov/2.0+column*(m_fov/(m_ncol-1)));
- };
- arma::vec SV(int row, int column, const arma::vec& attitude, double z) const{
- arma::vec rv(3);
- arma::Mat<double> rotM(3,3);
- rotM=getRz(attitude(2))*getRy(attitude(1))*getRx(attitude(0));
- rv=rotM*getIV(row,column);
- double height=rv(2);
- rv=rv*z/height;
- return rv;
- };
- arma::vec scanAngle_PB(int row, int column) const{
- arma::vec alpha(2);
- double r=corr_along(column);
- double theAx=(column<m_ppx) ? m_dx*0.000001*(m_ppx-(column+corr_across(column))) : m_dx*0.000001*(m_ppx-(column-corr_across(column)));
- double theAy=m_dy*0.000001*(m_ppy-r);
- alpha(0)=(getModel()==FRAME) ? atan(theAx/getF()) : atan(-theAx/getF());
- alpha(1)=atan(-theAy/getF());
- return alpha;
- };
- //geocentric coordinate system is right-handed, orthogogal Cartesian system with its origin at the centre of the earth. The direct Helmert transformation from lat/lon and ellipsoid height to geocentric x,y,z is returned in posX, posY, posZ
- //Euclidian distance between points on sphere
- //Euclidian distance between two vectors
- double getDist(const arma::vec& v1, const arma::vec& v2) const{arma::vec dv=v1;dv-=v2;return arma::norm(dv,2);};
- double lox1(double lat1, double lon1, double lat2, double lon2) const{
- double result=(M(lat2)-M(lat1))/cos(Az(lat1,lon1,lat2,lon2));
- return result;
- };
- double lox2(double lat, double lon1, double lon2) const{
- double eps2=getEccentricity();
- eps2*=eps2;
- double sin2=sin(lat);
- sin2*=sin2;
- return m_spatialref.GetSemiMajor()*(lon2-lon1)*cos(lat)/sqrt(1-eps2*sin2);
- };
- //length of the meridian from Equator to point
- double M(double lat) const{
- double eps2=getEccentricity();
- eps2*=eps2;
- double eps4=eps2*eps2;
- double eps6=eps4*eps2;
- double t1=(1-eps2/4.0-3*eps4/64.0-5*eps6/256.0)*lat;
- double t2=(3*eps2/8.0+3*eps4/32.0+45*eps6/1024.0)*sin(2*lat);
- double t3=(15*eps4/256.0+45*eps6/1024.0)*sin(4*lat);
- double t4=35*eps6/3072.0*sin(6*lat);
- double semiMajor=m_spatialref.GetSemiMajor();
- return semiMajor*(t1-t2+t3-t4);
- };
- double atan33(double x, double y) const{
- double angle=atan2(y,x);
- while(angle<0)
- angle+=2*PI;
- angle=PI/2.0+2*PI-angle;
- while(angle>2*PI)
- angle-=2*PI;
- return angle;
- };
- //The azimuth from North (clockwise) of the loxodrome from point 1 to point 2 will than be given by:
- double Az(double lat1, double lon1, double lat2, double lon2) const{
- double x1=lon1;
- double x2=lon2;
- double y1=lat1;
- double y2=lat2;
- getMxy(x1,y1);
- getMxy(x2,y2);
- return atan33(x2-x1,y2-y1);
- };
- double corr_along(double c) const{
- return 0;
- };
- double corr_across(double c) const{
- double result=0;
- if(m_polynome.size()){
- double absc=1;
- result+=absc*m_polynome[0];
- for(int degree=1;degree<m_polynome.size();++degree){
- absc*=sqrt((c-m_ppx)*(c-m_ppx));
- result+=absc*m_polynome[degree];
- }
- }
- return result;
- };
- arma::Mat<double> getRx(double theta) const{
- arma::Mat<double> Rx(3,3);
- Rx(0,0)=1;
- Rx(0,1)=0;
- Rx(0,2)=0;
- Rx(1,0)=0;
- Rx(1,1)=cos(theta);
- Rx(1,2)=-sin(theta);
- Rx(2,0)=0;
- Rx(2,1)=sin(theta);
- Rx(2,2)=cos(theta);
- return Rx;
- };
- arma::Mat<double> getRy(double theta) const{
- arma::Mat<double> Ry(3,3);
- Ry(0,0)=cos(theta);
- Ry(0,1)=0;
- Ry(0,2)=sin(theta);
- Ry(1,0)=0;
- Ry(1,1)=1;
- Ry(1,2)=0;
- Ry(2,0)=-sin(theta);
- Ry(2,1)=0;
- Ry(2,2)=cos(theta);
- return Ry;
- };
-
- arma::Mat<double> getRz(double theta) const{
- arma::Mat<double> Rz(3,3);
- Rz(0,0)=cos(theta);
- Rz(0,1)=-sin(theta);
- Rz(0,2)=0;
- Rz(1,0)=sin(theta);
- Rz(1,1)=cos(theta);
- Rz(1,2)=0;
- Rz(2,0)=0;
- Rz(2,1)=0;
- Rz(2,2)=1;
- return Rz;
- };
- int m_model;
- OGRSpatialReference m_spatialref;
- //Ae=m_spatialref.GetSemiMajor()
- //Be=m_spatialref.GetSemiMinor()
- //f_earth=1.0/m_spatialref.GetInvFlattening()
- std::string m_datum;
- double m_ppx;
- double m_ppy;
- double m_fc;
- double m_fov;
- int m_ncol;
- int m_nrow;
- double m_dx;
- double m_dy;
- arma::vec m_bcpos;
- arma::vec m_bcatt;
- std::vector<double> m_polynome;
- };
-}
-#endif //_SENSORMODEL_
diff --git a/src/models/dataSpec_P5B.f90 b/src/models/dataSpec_P5B.f90
deleted file mode 100644
index fbc4bd4..0000000
--- a/src/models/dataSpec_P5B.f90
+++ /dev/null
@@ -1,1139 +0,0 @@
-! ********************************************************************************
-! dataSpec_P5B.f90
-! ********************************************************************************
-! lambda = lambdalength (nm)
-! refractive = refractive index of leaf material
-! k_Cab = specific absorption coefficient of chlorophyll (a+b) (cm2.�g-1)
-! k_Car = specific absorption coefficient of carotenoids (cm2.�g-1)
-! k_Cw = specific absorption coefficient of water (cm-1)
-! k_Cm = specific absorption coefficient of dry matter (cm2.g-1)
-! k_Bp = specific absorption coefficient of brown pigments (arbitrary units)
-! ********************************************************************************
-! F�ret J.B, Fran�ois C, Asner G.P, Gitelson A.A, Martin R.E, Bidel L.P.R,
-! Ustin S.L, le Maire G, Jacquemoud S. (2008), PROSPECT-4 and 5: Advances in the
-! leaf optical properties model separating photosynthetic pigments, Remote Sensing
-! of Environment, 112:3030-3043.
-! The specific absorption coefficient of brown pigments is provided by F. Baret
-! (EMMAH, INRA Avignon, baret at avignon.inra.fr) and used with his autorization.
-! ********************************************************************************
-! version 5.02 (25 July 2011)
-! ********************************************************************************
-
-module dataSpec_P5B
-
-implicit none
-integer, parameter :: nw=2101
-integer i, lambda(nw)
-real(8) refractive(nw), k_Cab(nw), k_Car(nw), k_Brown(nw), k_Cw(nw), k_Cm(nw)
-
-
-! ********************************************************************************
-! Wavelength
-! ********************************************************************************
- data (lambda(i),i=1,100)/&
-400,401,402,403,404,405,406,407,408,409,&
-410,411,412,413,414,415,416,417,418,419,&
-420,421,422,423,424,425,426,427,428,429,&
-430,431,432,433,434,435,436,437,438,439,&
-440,441,442,443,444,445,446,447,448,449,&
-450,451,452,453,454,455,456,457,458,459,&
-460,461,462,463,464,465,466,467,468,469,&
-470,471,472,473,474,475,476,477,478,479,&
-480,481,482,483,484,485,486,487,488,489,&
-490,491,492,493,494,495,496,497,498,499/
- data (lambda(i),i=101,200)/&
-500,501,502,503,504,505,506,507,508,509,&
-510,511,512,513,514,515,516,517,518,519,&
-520,521,522,523,524,525,526,527,528,529,&
-530,531,532,533,534,535,536,537,538,539,&
-540,541,542,543,544,545,546,547,548,549,&
-550,551,552,553,554,555,556,557,558,559,&
-560,561,562,563,564,565,566,567,568,569,&
-570,571,572,573,574,575,576,577,578,579,&
-580,581,582,583,584,585,586,587,588,589,&
-590,591,592,593,594,595,596,597,598,599/
- data (lambda(i),i=201,300)/&
-600,601,602,603,604,605,606,607,608,609,&
-610,611,612,613,614,615,616,617,618,619,&
-620,621,622,623,624,625,626,627,628,629,&
-630,631,632,633,634,635,636,637,638,639,&
-640,641,642,643,644,645,646,647,648,649,&
-650,651,652,653,654,655,656,657,658,659,&
-660,661,662,663,664,665,666,667,668,669,&
-670,671,672,673,674,675,676,677,678,679,&
-680,681,682,683,684,685,686,687,688,689,&
-690,691,692,693,694,695,696,697,698,699/
- data (lambda(i),i=301,400)/&
-700,701,702,703,704,705,706,707,708,709,&
-710,711,712,713,714,715,716,717,718,719,&
-720,721,722,723,724,725,726,727,728,729,&
-730,731,732,733,734,735,736,737,738,739,&
-740,741,742,743,744,745,746,747,748,749,&
-750,751,752,753,754,755,756,757,758,759,&
-760,761,762,763,764,765,766,767,768,769,&
-770,771,772,773,774,775,776,777,778,779,&
-780,781,782,783,784,785,786,787,788,789,&
-790,791,792,793,794,795,796,797,798,799/
- data (lambda(i),i=401,500)/&
-800,801,802,803,804,805,806,807,808,809,&
-810,811,812,813,814,815,816,817,818,819,&
-820,821,822,823,824,825,826,827,828,829,&
-830,831,832,833,834,835,836,837,838,839,&
-840,841,842,843,844,845,846,847,848,849,&
-850,851,852,853,854,855,856,857,858,859,&
-860,861,862,863,864,865,866,867,868,869,&
-870,871,872,873,874,875,876,877,878,879,&
-880,881,882,883,884,885,886,887,888,889,&
-890,891,892,893,894,895,896,897,898,899/
- data (lambda(i),i=501,600)/&
-900,901,902,903,904,905,906,907,908,909,&
-910,911,912,913,914,915,916,917,918,919,&
-920,921,922,923,924,925,926,927,928,929,&
-930,931,932,933,934,935,936,937,938,939,&
-940,941,942,943,944,945,946,947,948,949,&
-950,951,952,953,954,955,956,957,958,959,&
-960,961,962,963,964,965,966,967,968,969,&
-970,971,972,973,974,975,976,977,978,979,&
-980,981,982,983,984,985,986,987,988,989,&
-990,991,992,993,994,995,996,997,998,999/
- data (lambda(i),i=601,700)/&
-1000,1001,1002,1003,1004,1005,1006,1007,1008,1009,&
-1010,1011,1012,1013,1014,1015,1016,1017,1018,1019,&
-1020,1021,1022,1023,1024,1025,1026,1027,1028,1029,&
-1030,1031,1032,1033,1034,1035,1036,1037,1038,1039,&
-1040,1041,1042,1043,1044,1045,1046,1047,1048,1049,&
-1050,1051,1052,1053,1054,1055,1056,1057,1058,1059,&
-1060,1061,1062,1063,1064,1065,1066,1067,1068,1069,&
-1070,1071,1072,1073,1074,1075,1076,1077,1078,1079,&
-1080,1081,1082,1083,1084,1085,1086,1087,1088,1089,&
-1090,1091,1092,1093,1094,1095,1096,1097,1098,1099/
- data (lambda(i),i=701,800)/&
-1100,1101,1102,1103,1104,1105,1106,1107,1108,1109,&
-1110,1111,1112,1113,1114,1115,1116,1117,1118,1119,&
-1120,1121,1122,1123,1124,1125,1126,1127,1128,1129,&
-1130,1131,1132,1133,1134,1135,1136,1137,1138,1139,&
-1140,1141,1142,1143,1144,1145,1146,1147,1148,1149,&
-1150,1151,1152,1153,1154,1155,1156,1157,1158,1159,&
-1160,1161,1162,1163,1164,1165,1166,1167,1168,1169,&
-1170,1171,1172,1173,1174,1175,1176,1177,1178,1179,&
-1180,1181,1182,1183,1184,1185,1186,1187,1188,1189,&
-1190,1191,1192,1193,1194,1195,1196,1197,1198,1199/
- data (lambda(i),i=801,900)/&
-1200,1201,1202,1203,1204,1205,1206,1207,1208,1209,&
-1210,1211,1212,1213,1214,1215,1216,1217,1218,1219,&
-1220,1221,1222,1223,1224,1225,1226,1227,1228,1229,&
-1230,1231,1232,1233,1234,1235,1236,1237,1238,1239,&
-1240,1241,1242,1243,1244,1245,1246,1247,1248,1249,&
-1250,1251,1252,1253,1254,1255,1256,1257,1258,1259,&
-1260,1261,1262,1263,1264,1265,1266,1267,1268,1269,&
-1270,1271,1272,1273,1274,1275,1276,1277,1278,1279,&
-1280,1281,1282,1283,1284,1285,1286,1287,1288,1289,&
-1290,1291,1292,1293,1294,1295,1296,1297,1298,1299/
- data (lambda(i),i=901,1000)/&
-1300,1301,1302,1303,1304,1305,1306,1307,1308,1309,&
-1310,1311,1312,1313,1314,1315,1316,1317,1318,1319,&
-1320,1321,1322,1323,1324,1325,1326,1327,1328,1329,&
-1330,1331,1332,1333,1334,1335,1336,1337,1338,1339,&
-1340,1341,1342,1343,1344,1345,1346,1347,1348,1349,&
-1350,1351,1352,1353,1354,1355,1356,1357,1358,1359,&
-1360,1361,1362,1363,1364,1365,1366,1367,1368,1369,&
-1370,1371,1372,1373,1374,1375,1376,1377,1378,1379,&
-1380,1381,1382,1383,1384,1385,1386,1387,1388,1389,&
-1390,1391,1392,1393,1394,1395,1396,1397,1398,1399/
- data (lambda(i),i=1001,1100)/&
-1400,1401,1402,1403,1404,1405,1406,1407,1408,1409,&
-1410,1411,1412,1413,1414,1415,1416,1417,1418,1419,&
-1420,1421,1422,1423,1424,1425,1426,1427,1428,1429,&
-1430,1431,1432,1433,1434,1435,1436,1437,1438,1439,&
-1440,1441,1442,1443,1444,1445,1446,1447,1448,1449,&
-1450,1451,1452,1453,1454,1455,1456,1457,1458,1459,&
-1460,1461,1462,1463,1464,1465,1466,1467,1468,1469,&
-1470,1471,1472,1473,1474,1475,1476,1477,1478,1479,&
-1480,1481,1482,1483,1484,1485,1486,1487,1488,1489,&
-1490,1491,1492,1493,1494,1495,1496,1497,1498,1499/
- data (lambda(i),i=1101,1200)/&
-1500,1501,1502,1503,1504,1505,1506,1507,1508,1509,&
-1510,1511,1512,1513,1514,1515,1516,1517,1518,1519,&
-1520,1521,1522,1523,1524,1525,1526,1527,1528,1529,&
-1530,1531,1532,1533,1534,1535,1536,1537,1538,1539,&
-1540,1541,1542,1543,1544,1545,1546,1547,1548,1549,&
-1550,1551,1552,1553,1554,1555,1556,1557,1558,1559,&
-1560,1561,1562,1563,1564,1565,1566,1567,1568,1569,&
-1570,1571,1572,1573,1574,1575,1576,1577,1578,1579,&
-1580,1581,1582,1583,1584,1585,1586,1587,1588,1589,&
-1590,1591,1592,1593,1594,1595,1596,1597,1598,1599/
- data (lambda(i),i=1201,1300)/&
-1600,1601,1602,1603,1604,1605,1606,1607,1608,1609,&
-1610,1611,1612,1613,1614,1615,1616,1617,1618,1619,&
-1620,1621,1622,1623,1624,1625,1626,1627,1628,1629,&
-1630,1631,1632,1633,1634,1635,1636,1637,1638,1639,&
-1640,1641,1642,1643,1644,1645,1646,1647,1648,1649,&
-1650,1651,1652,1653,1654,1655,1656,1657,1658,1659,&
-1660,1661,1662,1663,1664,1665,1666,1667,1668,1669,&
-1670,1671,1672,1673,1674,1675,1676,1677,1678,1679,&
-1680,1681,1682,1683,1684,1685,1686,1687,1688,1689,&
-1690,1691,1692,1693,1694,1695,1696,1697,1698,1699/
- data (lambda(i),i=1301,1400)/&
-1700,1701,1702,1703,1704,1705,1706,1707,1708,1709,&
-1710,1711,1712,1713,1714,1715,1716,1717,1718,1719,&
-1720,1721,1722,1723,1724,1725,1726,1727,1728,1729,&
-1730,1731,1732,1733,1734,1735,1736,1737,1738,1739,&
-1740,1741,1742,1743,1744,1745,1746,1747,1748,1749,&
-1750,1751,1752,1753,1754,1755,1756,1757,1758,1759,&
-1760,1761,1762,1763,1764,1765,1766,1767,1768,1769,&
-1770,1771,1772,1773,1774,1775,1776,1777,1778,1779,&
-1780,1781,1782,1783,1784,1785,1786,1787,1788,1789,&
-1790,1791,1792,1793,1794,1795,1796,1797,1798,1799/
- data (lambda(i),i=1401,1500)/&
-1800,1801,1802,1803,1804,1805,1806,1807,1808,1809,&
-1810,1811,1812,1813,1814,1815,1816,1817,1818,1819,&
-1820,1821,1822,1823,1824,1825,1826,1827,1828,1829,&
-1830,1831,1832,1833,1834,1835,1836,1837,1838,1839,&
-1840,1841,1842,1843,1844,1845,1846,1847,1848,1849,&
-1850,1851,1852,1853,1854,1855,1856,1857,1858,1859,&
-1860,1861,1862,1863,1864,1865,1866,1867,1868,1869,&
-1870,1871,1872,1873,1874,1875,1876,1877,1878,1879,&
-1880,1881,1882,1883,1884,1885,1886,1887,1888,1889,&
-1890,1891,1892,1893,1894,1895,1896,1897,1898,1899/
- data (lambda(i),i=1501,1600)/&
-1900,1901,1902,1903,1904,1905,1906,1907,1908,1909,&
-1910,1911,1912,1913,1914,1915,1916,1917,1918,1919,&
-1920,1921,1922,1923,1924,1925,1926,1927,1928,1929,&
-1930,1931,1932,1933,1934,1935,1936,1937,1938,1939,&
-1940,1941,1942,1943,1944,1945,1946,1947,1948,1949,&
-1950,1951,1952,1953,1954,1955,1956,1957,1958,1959,&
-1960,1961,1962,1963,1964,1965,1966,1967,1968,1969,&
-1970,1971,1972,1973,1974,1975,1976,1977,1978,1979,&
-1980,1981,1982,1983,1984,1985,1986,1987,1988,1989,&
-1990,1991,1992,1993,1994,1995,1996,1997,1998,1999/
- data (lambda(i),i=1601,1700)/&
-2000,2001,2002,2003,2004,2005,2006,2007,2008,2009,&
-2010,2011,2012,2013,2014,2015,2016,2017,2018,2019,&
-2020,2021,2022,2023,2024,2025,2026,2027,2028,2029,&
-2030,2031,2032,2033,2034,2035,2036,2037,2038,2039,&
-2040,2041,2042,2043,2044,2045,2046,2047,2048,2049,&
-2050,2051,2052,2053,2054,2055,2056,2057,2058,2059,&
-2060,2061,2062,2063,2064,2065,2066,2067,2068,2069,&
-2070,2071,2072,2073,2074,2075,2076,2077,2078,2079,&
-2080,2081,2082,2083,2084,2085,2086,2087,2088,2089,&
-2090,2091,2092,2093,2094,2095,2096,2097,2098,2099/
- data (lambda(i),i=1701,1800)/&
-2100,2101,2102,2103,2104,2105,2106,2107,2108,2109,&
-2110,2111,2112,2113,2114,2115,2116,2117,2118,2119,&
-2120,2121,2122,2123,2124,2125,2126,2127,2128,2129,&
-2130,2131,2132,2133,2134,2135,2136,2137,2138,2139,&
-2140,2141,2142,2143,2144,2145,2146,2147,2148,2149,&
-2150,2151,2152,2153,2154,2155,2156,2157,2158,2159,&
-2160,2161,2162,2163,2164,2165,2166,2167,2168,2169,&
-2170,2171,2172,2173,2174,2175,2176,2177,2178,2179,&
-2180,2181,2182,2183,2184,2185,2186,2187,2188,2189,&
-2190,2191,2192,2193,2194,2195,2196,2197,2198,2199/
- data (lambda(i),i=1801,1900)/&
-2200,2201,2202,2203,2204,2205,2206,2207,2208,2209,&
-2210,2211,2212,2213,2214,2215,2216,2217,2218,2219,&
-2220,2221,2222,2223,2224,2225,2226,2227,2228,2229,&
-2230,2231,2232,2233,2234,2235,2236,2237,2238,2239,&
-2240,2241,2242,2243,2244,2245,2246,2247,2248,2249,&
-2250,2251,2252,2253,2254,2255,2256,2257,2258,2259,&
-2260,2261,2262,2263,2264,2265,2266,2267,2268,2269,&
-2270,2271,2272,2273,2274,2275,2276,2277,2278,2279,&
-2280,2281,2282,2283,2284,2285,2286,2287,2288,2289,&
-2290,2291,2292,2293,2294,2295,2296,2297,2298,2299/
- data (lambda(i),i=1901,2000)/&
-2300,2301,2302,2303,2304,2305,2306,2307,2308,2309,&
-2310,2311,2312,2313,2314,2315,2316,2317,2318,2319,&
-2320,2321,2322,2323,2324,2325,2326,2327,2328,2329,&
-2330,2331,2332,2333,2334,2335,2336,2337,2338,2339,&
-2340,2341,2342,2343,2344,2345,2346,2347,2348,2349,&
-2350,2351,2352,2353,2354,2355,2356,2357,2358,2359,&
-2360,2361,2362,2363,2364,2365,2366,2367,2368,2369,&
-2370,2371,2372,2373,2374,2375,2376,2377,2378,2379,&
-2380,2381,2382,2383,2384,2385,2386,2387,2388,2389,&
-2390,2391,2392,2393,2394,2395,2396,2397,2398,2399/
- data (lambda(i),i=2001,2101)/&
-2400,2401,2402,2403,2404,2405,2406,2407,2408,2409,&
-2410,2411,2412,2413,2414,2415,2416,2417,2418,2419,&
-2420,2421,2422,2423,2424,2425,2426,2427,2428,2429,&
-2430,2431,2432,2433,2434,2435,2436,2437,2438,2439,&
-2440,2441,2442,2443,2444,2445,2446,2447,2448,2449,&
-2450,2451,2452,2453,2454,2455,2456,2457,2458,2459,&
-2460,2461,2462,2463,2464,2465,2466,2467,2468,2469,&
-2470,2471,2472,2473,2474,2475,2476,2477,2478,2479,&
-2480,2481,2482,2483,2484,2485,2486,2487,2488,2489,&
-2490,2491,2492,2493,2494,2495,2496,2497,2498,2499,&
-2500./
-
-! ********************************************************************************
-! Refractive index
-! ********************************************************************************
-
- data (refractive(i),i=1,100)/&
-1.4955,1.4958,1.4960,1.4964,1.4971,1.4978,1.4986,1.4995,1.5006,1.5016,&
-1.5024,1.5031,1.5041,1.5052,1.5065,1.5076,1.5088,1.5099,1.5110,1.5122,&
-1.5134,1.5144,1.5156,1.5168,1.5181,1.5195,1.5206,1.5212,1.5219,1.5228,&
-1.5235,1.5241,1.5249,1.5256,1.5262,1.5267,1.5272,1.5276,1.5279,1.5282,&
-1.5286,1.5290,1.5292,1.5294,1.5294,1.5295,1.5293,1.5291,1.5289,1.5286,&
-1.5282,1.5278,1.5273,1.5268,1.5262,1.5256,1.5249,1.5243,1.5237,1.5230,&
-1.5224,1.5217,1.5210,1.5203,1.5196,1.5189,1.5184,1.5178,1.5173,1.5170,&
-1.5166,1.5163,1.5159,1.5156,1.5153,1.5150,1.5148,1.5145,1.5142,1.5139,&
-1.5135,1.5131,1.5125,1.5119,1.5111,1.5103,1.5092,1.5080,1.5067,1.5051,&
-1.5034,1.5015,1.4994,1.4971,1.4947,1.4921,1.4894,1.4866,1.4837,1.4807/
- data (refractive(i),i=101,200)/&
-1.4776,1.4744,1.4712,1.4680,1.4647,1.4615,1.4583,1.4553,1.4523,1.4495,&
-1.4468,1.4443,1.4419,1.4397,1.4377,1.4360,1.4344,1.4331,1.4319,1.4310,&
-1.4302,1.4296,1.4291,1.4288,1.4286,1.4286,1.4286,1.4287,1.4288,1.4290,&
-1.4292,1.4294,1.4296,1.4298,1.4300,1.4302,1.4303,1.4305,1.4306,1.4307,&
-1.4308,1.4309,1.4309,1.4310,1.4310,1.4310,1.4309,1.4308,1.4307,1.4305,&
-1.4303,1.4301,1.4298,1.4294,1.4290,1.4285,1.4280,1.4274,1.4268,1.4262,&
-1.4255,1.4248,1.4241,1.4234,1.4227,1.4220,1.4213,1.4206,1.4199,1.4192,&
-1.4185,1.4179,1.4173,1.4167,1.4161,1.4156,1.4151,1.4146,1.4141,1.4137,&
-1.4133,1.4129,1.4126,1.4123,1.4120,1.4117,1.4114,1.4112,1.4109,1.4107,&
-1.4105,1.4103,1.4101,1.4099,1.4097,1.4095,1.4093,1.4091,1.4089,1.4087/
- data (refractive(i),i=201,300)/&
-1.4086,1.4084,1.4082,1.4081,1.4079,1.4078,1.4077,1.4076,1.4075,1.4075,&
-1.4075,1.4075,1.4075,1.4076,1.4077,1.4077,1.4078,1.4079,1.4080,1.4081,&
-1.4082,1.4083,1.4083,1.4083,1.4083,1.4082,1.4082,1.4081,1.4079,1.4078,&
-1.4077,1.4076,1.4076,1.4076,1.4078,1.4080,1.4083,1.4087,1.4093,1.4100,&
-1.4108,1.4118,1.4128,1.4140,1.4153,1.4166,1.4179,1.4193,1.4206,1.4218,&
-1.4229,1.4240,1.4250,1.4261,1.4272,1.4284,1.4298,1.4315,1.4335,1.4357,&
-1.4383,1.4411,1.4442,1.4475,1.4509,1.4544,1.4579,1.4614,1.4649,1.4682,&
-1.4714,1.4744,1.4773,1.4802,1.4830,1.4857,1.4884,1.4909,1.4932,1.4951,&
-1.4963,1.4967,1.4958,1.4936,1.4899,1.4847,1.4780,1.4702,1.4618,1.4530,&
-1.4445,1.4366,1.4296,1.4237,1.4190,1.4154,1.4129,1.4112,1.4104,1.4102/
- data (refractive(i),i=301,400)/&
-1.4106,1.4113,1.4124,1.4137,1.4152,1.4168,1.4186,1.4204,1.4222,1.4242,&
-1.4261,1.4281,1.4301,1.4321,1.4341,1.4361,1.4382,1.4402,1.4422,1.4442,&
-1.4462,1.4481,1.4500,1.4518,1.4536,1.4554,1.4571,1.4587,1.4603,1.4618,&
-1.4632,1.4646,1.4659,1.4671,1.4683,1.4694,1.4704,1.4713,1.4721,1.4729,&
-1.4736,1.4742,1.4747,1.4752,1.4756,1.4759,1.4762,1.4764,1.4765,1.4766,&
-1.4766,1.4767,1.4767,1.4767,1.4767,1.4766,1.4766,1.4765,1.4764,1.4763,&
-1.4762,1.4761,1.4759,1.4758,1.4756,1.4754,1.4752,1.4751,1.4748,1.4746,&
-1.4744,1.4742,1.4739,1.4737,1.4734,1.4732,1.4729,1.4727,1.4725,1.4722,&
-1.4720,1.4718,1.4716,1.4714,1.4712,1.4710,1.4708,1.4706,1.4703,1.4701,&
-1.4700,1.4698,1.4696,1.4694,1.4692,1.4690,1.4688,1.4686,1.4684,1.4682/
- data (refractive(i),i=401,500)/&
-1.4681,1.4679,1.4677,1.4675,1.4674,1.4672,1.4670,1.4668,1.4666,1.4665,&
-1.4663,1.4662,1.4660,1.4659,1.4657,1.4655,1.4654,1.4652,1.4651,1.4649,&
-1.4648,1.4646,1.4645,1.4643,1.4641,1.4639,1.4638,1.4636,1.4635,1.4633,&
-1.4631,1.4630,1.4628,1.4627,1.4626,1.4625,1.4623,1.4622,1.4620,1.4619,&
-1.4617,1.4616,1.4614,1.4613,1.4611,1.4610,1.4608,1.4607,1.4605,1.4604,&
-1.4602,1.4601,1.4600,1.4599,1.4597,1.4596,1.4595,1.4594,1.4592,1.4591,&
-1.4589,1.4587,1.4586,1.4584,1.4583,1.4582,1.4582,1.4581,1.4580,1.4579,&
-1.4577,1.4576,1.4574,1.4572,1.4571,1.4569,1.4568,1.4566,1.4565,1.4563,&
-1.4562,1.4560,1.4559,1.4558,1.4556,1.4554,1.4553,1.4551,1.4550,1.4548,&
-1.4547,1.4545,1.4544,1.4542,1.4541,1.4539,1.4538,1.4536,1.4535,1.4533/
- data (refractive(i),i=501,600)/&
-1.4532,1.4531,1.4529,1.4528,1.4526,1.4525,1.4524,1.4522,1.4521,1.4519,&
-1.4518,1.4516,1.4515,1.4513,1.4512,1.4511,1.4509,1.4508,1.4506,1.4505,&
-1.4503,1.4502,1.4500,1.4498,1.4497,1.4495,1.4493,1.4492,1.4491,1.4489,&
-1.4488,1.4487,1.4485,1.4483,1.4481,1.4479,1.4478,1.4477,1.4476,1.4474,&
-1.4472,1.4469,1.4466,1.4465,1.4463,1.4462,1.4461,1.4461,1.4461,1.4460,&
-1.4457,1.4455,1.4453,1.4451,1.4450,1.4449,1.4448,1.4447,1.4445,1.4445,&
-1.4445,1.4442,1.4440,1.4438,1.4437,1.4436,1.4436,1.4435,1.4434,1.4434,&
-1.4434,1.4434,1.4432,1.4430,1.4428,1.4427,1.4426,1.4425,1.4423,1.4421,&
-1.4420,1.4418,1.4416,1.4415,1.4414,1.4414,1.4414,1.4415,1.4416,1.4416,&
-1.4415,1.4414,1.4413,1.4411,1.4409,1.4407,1.4405,1.4403,1.4401,1.4400/
- data (refractive(i),i=601,700)/&
-1.4399,1.4397,1.4395,1.4394,1.4393,1.4392,1.4392,1.4391,1.4391,1.4390,&
-1.4389,1.4387,1.4386,1.4385,1.4384,1.4384,1.4383,1.4381,1.4379,1.4378,&
-1.4376,1.4374,1.4372,1.4371,1.4370,1.4369,1.4368,1.4367,1.4365,1.4365,&
-1.4364,1.4364,1.4363,1.4363,1.4363,1.4363,1.4363,1.4362,1.4361,1.4360,&
-1.4359,1.4357,1.4356,1.4354,1.4352,1.4351,1.4351,1.4350,1.4349,1.4349,&
-1.4348,1.4347,1.4345,1.4343,1.4342,1.4342,1.4342,1.4341,1.4340,1.4339,&
-1.4338,1.4337,1.4335,1.4334,1.4333,1.4333,1.4334,1.4333,1.4333,1.4332,&
-1.4330,1.4329,1.4328,1.4327,1.4325,1.4324,1.4322,1.4321,1.4319,1.4318,&
-1.4317,1.4316,1.4315,1.4314,1.4314,1.4313,1.4312,1.4310,1.4309,1.4307,&
-1.4305,1.4304,1.4303,1.4302,1.4302,1.4301,1.4300,1.4297,1.4295,1.4294/
- data (refractive(i),i=701,800)/&
-1.4293,1.4292,1.4291,1.4290,1.4289,1.4288,1.4286,1.4284,1.4283,1.4283,&
-1.4284,1.4284,1.4284,1.4284,1.4283,1.4282,1.4282,1.4280,1.4278,1.4276,&
-1.4275,1.4274,1.4272,1.4270,1.4268,1.4267,1.4267,1.4266,1.4265,1.4264,&
-1.4263,1.4262,1.4260,1.4259,1.4258,1.4257,1.4256,1.4256,1.4255,1.4254,&
-1.4252,1.4250,1.4249,1.4248,1.4247,1.4244,1.4241,1.4239,1.4237,1.4235,&
-1.4234,1.4233,1.4231,1.4230,1.4228,1.4226,1.4224,1.4224,1.4223,1.4223,&
-1.4222,1.4222,1.4221,1.4221,1.4221,1.4219,1.4217,1.4216,1.4215,1.4215,&
-1.4214,1.4212,1.4210,1.4208,1.4207,1.4206,1.4206,1.4206,1.4206,1.4206,&
-1.4205,1.4203,1.4200,1.4199,1.4198,1.4199,1.4201,1.4204,1.4205,1.4205,&
-1.4203,1.4200,1.4196,1.4192,1.4189,1.4187,1.4187,1.4188,1.4188,1.4187/
- data (refractive(i),i=801,900)/&
-1.4186,1.4185,1.4182,1.4181,1.4179,1.4178,1.4177,1.4177,1.4176,1.4175,&
-1.4173,1.4173,1.4172,1.4171,1.4170,1.4168,1.4167,1.4165,1.4164,1.4163,&
-1.4161,1.4158,1.4156,1.4155,1.4154,1.4153,1.4151,1.4151,1.4151,1.4150,&
-1.4148,1.4147,1.4146,1.4145,1.4144,1.4143,1.4141,1.4139,1.4136,1.4134,&
-1.4133,1.4130,1.4128,1.4128,1.4128,1.4127,1.4126,1.4124,1.4122,1.4122,&
-1.4121,1.4121,1.4121,1.4119,1.4118,1.4116,1.4115,1.4113,1.4112,1.4112,&
-1.4111,1.4111,1.4110,1.4108,1.4106,1.4104,1.4102,1.4102,1.4102,1.4103,&
-1.4103,1.4101,1.4099,1.4097,1.4095,1.4093,1.4091,1.4091,1.4091,1.4090,&
-1.4089,1.4088,1.4085,1.4083,1.4081,1.4080,1.4079,1.4078,1.4077,1.4078,&
-1.4078,1.4078,1.4077,1.4076,1.4075,1.4073,1.4070,1.4068,1.4066,1.4064/
- data (refractive(i),i=901,1000)/&
-1.4064,1.4062,1.4061,1.4059,1.4056,1.4055,1.4055,1.4055,1.4055,1.4056,&
-1.4057,1.4058,1.4058,1.4056,1.4055,1.4053,1.4051,1.4050,1.4049,1.4048,&
-1.4046,1.4045,1.4044,1.4044,1.4044,1.4045,1.4045,1.4043,1.4042,1.4040,&
-1.4040,1.4041,1.4041,1.4040,1.4036,1.4034,1.4032,1.4031,1.4030,1.4031,&
-1.4033,1.4033,1.4033,1.4031,1.4028,1.4023,1.4017,1.4013,1.4010,1.4010,&
-1.4008,1.4007,1.4004,1.4002,1.4000,1.4002,1.4004,1.4009,1.4011,1.4015,&
-1.4016,1.4015,1.4010,1.4002,1.3996,1.3992,1.3989,1.3992,1.3994,1.3997,&
-1.4000,1.4000,1.4000,1.4003,1.4002,1.3999,1.4002,1.4004,1.4008,1.4009,&
-1.4019,1.4015,1.4010,1.4000,1.4001,1.3988,1.3978,1.3971,1.3960,1.3967,&
-1.3958,1.3955,1.3952,1.3953,1.3945,1.3960,1.3961,1.3950,1.3946,1.3942/
- data (refractive(i),i=1001,1100)/&
-1.3938,1.3935,1.3926,1.3926,1.3930,1.3943,1.3943,1.3944,1.3944,1.3937,&
-1.3933,1.3934,1.3934,1.3935,1.3933,1.3931,1.3928,1.3925,1.3920,1.3918,&
-1.3915,1.3910,1.3904,1.3898,1.3892,1.3880,1.3872,1.3865,1.3858,1.3852,&
-1.3843,1.3838,1.3829,1.3823,1.3810,1.3801,1.3794,1.3785,1.3776,1.3771,&
-1.3765,1.3758,1.3751,1.3743,1.3733,1.3725,1.3717,1.3711,1.3707,1.3703,&
-1.3700,1.3700,1.3697,1.3692,1.3685,1.3680,1.3675,1.3670,1.3667,1.3662,&
-1.3657,1.3659,1.3658,1.3656,1.3656,1.3654,1.3650,1.3653,1.3651,1.3649,&
-1.3650,1.3646,1.3641,1.3644,1.3640,1.3644,1.3643,1.3644,1.3647,1.3646,&
-1.3649,1.3644,1.3641,1.3635,1.3635,1.3633,1.3631,1.3635,1.3635,1.3643,&
-1.3642,1.3646,1.3643,1.3642,1.3637,1.3640,1.3638,1.3640,1.3640,1.3640/
- data (refractive(i),i=1101,1200)/&
-1.3646,1.3645,1.3650,1.3650,1.3652,1.3657,1.3656,1.3660,1.3657,1.3661,&
-1.3665,1.3667,1.3670,1.3666,1.3669,1.3673,1.3673,1.3680,1.3686,1.3686,&
-1.3689,1.3693,1.3690,1.3694,1.3697,1.3695,1.3697,1.3698,1.3696,1.3698,&
-1.3700,1.3706,1.3705,1.3709,1.3712,1.3717,1.3718,1.3721,1.3721,1.3724,&
-1.3723,1.3721,1.3723,1.3723,1.3725,1.3723,1.3723,1.3727,1.3729,1.3730,&
-1.3730,1.3730,1.3731,1.3734,1.3734,1.3736,1.3737,1.3738,1.3735,1.3735,&
-1.3734,1.3735,1.3736,1.3738,1.3737,1.3738,1.3736,1.3737,1.3735,1.3735,&
-1.3735,1.3734,1.3734,1.3734,1.3735,1.3732,1.3733,1.3732,1.3732,1.3735,&
-1.3736,1.3739,1.3742,1.3741,1.3738,1.3737,1.3737,1.3736,1.3735,1.3737,&
-1.3738,1.3741,1.3743,1.3744,1.3744,1.3743,1.3740,1.3740,1.3741,1.3741/
- data (refractive(i),i=1201,1300)/&
-1.3742,1.3742,1.3741,1.3743,1.3743,1.3743,1.3742,1.3741,1.3740,1.3740,&
-1.3738,1.3736,1.3734,1.3733,1.3731,1.3730,1.3730,1.3731,1.3730,1.3732,&
-1.3731,1.3730,1.3729,1.3728,1.3727,1.3726,1.3725,1.3726,1.3725,1.3726,&
-1.3726,1.3724,1.3719,1.3718,1.3716,1.3715,1.3712,1.3713,1.3711,1.3710,&
-1.3707,1.3708,1.3707,1.3708,1.3708,1.3710,1.3712,1.3715,1.3715,1.3716,&
-1.3714,1.3714,1.3712,1.3713,1.3709,1.3708,1.3704,1.3704,1.3700,1.3698,&
-1.3695,1.3694,1.3692,1.3690,1.3687,1.3687,1.3684,1.3685,1.3684,1.3683,&
-1.3683,1.3687,1.3681,1.3678,1.3674,1.3670,1.3665,1.3665,1.3664,1.3663,&
-1.3659,1.3659,1.3656,1.3655,1.3654,1.3655,1.3654,1.3654,1.3651,1.3651,&
-1.3648,1.3648,1.3647,1.3648,1.3647,1.3647,1.3646,1.3647,1.3646,1.3644/
- data (refractive(i),i=1301,1400)/&
-1.3640,1.3638,1.3635,1.3634,1.3632,1.3631,1.3628,1.3629,1.3627,1.3629,&
-1.3629,1.3627,1.3625,1.3625,1.3623,1.3622,1.3619,1.3620,1.3619,1.3619,&
-1.3616,1.3615,1.3612,1.3611,1.3611,1.3612,1.3613,1.3617,1.3617,1.3618,&
-1.3618,1.3617,1.3615,1.3616,1.3615,1.3616,1.3618,1.3617,1.3616,1.3616,&
-1.3613,1.3612,1.3612,1.3616,1.3616,1.3617,1.3618,1.3621,1.3623,1.3625,&
-1.3626,1.3626,1.3625,1.3627,1.3633,1.3637,1.3637,1.3636,1.3635,1.3633,&
-1.3630,1.3629,1.3632,1.3636,1.3641,1.3646,1.3651,1.3653,1.3648,1.3640,&
-1.3630,1.3620,1.3615,1.3608,1.3608,1.3611,1.3617,1.3628,1.3638,1.3645,&
-1.3649,1.3646,1.3637,1.3625,1.3613,1.3601,1.3595,1.3593,1.3599,1.3607,&
-1.3614,1.3623,1.3627,1.3627,1.3624,1.3622,1.3619,1.3616,1.3616,1.3617/
- data (refractive(i),i=1401,1500)/&
-1.3619,1.3621,1.3622,1.3626,1.3630,1.3635,1.3637,1.3639,1.3636,1.3632,&
-1.3629,1.3625,1.3623,1.3625,1.3628,1.3634,1.3638,1.3645,1.3649,1.3652,&
-1.3652,1.3650,1.3645,1.3642,1.3639,1.3638,1.3636,1.3634,1.3631,1.3632,&
-1.3631,1.3630,1.3624,1.3623,1.3623,1.3623,1.3626,1.3629,1.3629,1.3630,&
-1.3628,1.3622,1.3615,1.3607,1.3601,1.3597,1.3597,1.3599,1.3597,1.3597,&
-1.3595,1.3599,1.3593,1.3597,1.3602,1.3613,1.3625,1.3633,1.3638,1.3645,&
-1.3643,1.3648,1.3629,1.3632,1.3635,1.3645,1.3647,1.3662,1.3668,1.3684,&
-1.3653,1.3658,1.3663,1.3663,1.3671,1.3668,1.3680,1.3626,1.3620,1.3624,&
-1.3627,1.3617,1.3628,1.3562,1.3557,1.3558,1.3559,1.3555,1.3508,1.3504,&
-1.3511,1.3520,1.3525,1.3493,1.3503,1.3524,1.3546,1.3578,1.3556,1.3583/
- data (refractive(i),i=1501,1600)/&
-1.3601,1.3610,1.3595,1.3599,1.3606,1.3611,1.3598,1.3600,1.3601,1.3606,&
-1.3599,1.3600,1.3603,1.3608,1.3606,1.3613,1.3615,1.3614,1.3606,1.3600,&
-1.3591,1.3582,1.3573,1.3564,1.3556,1.3547,1.3538,1.3530,1.3522,1.3513,&
-1.3505,1.3497,1.3489,1.3481,1.3473,1.3465,1.3457,1.3450,1.3442,1.3434,&
-1.3427,1.3420,1.3412,1.3405,1.3398,1.3391,1.3384,1.3377,1.3371,1.3364,&
-1.3357,1.3351,1.3344,1.3338,1.3332,1.3326,1.3320,1.3314,1.3308,1.3302,&
-1.3297,1.3291,1.3286,1.3280,1.3275,1.3270,1.3265,1.3260,1.3255,1.3251,&
-1.3246,1.3242,1.3237,1.3233,1.3229,1.3225,1.3221,1.3217,1.3213,1.3209,&
-1.3206,1.3203,1.3199,1.3196,1.3193,1.3190,1.3187,1.3185,1.3182,1.3180,&
-1.3177,1.3175,1.3173,1.3171,1.3169,1.3168,1.3166,1.3165,1.3163,1.3162/
- data (refractive(i),i=1601,1700)/&
-1.3161,1.3160,1.3159,1.3159,1.3158,1.3157,1.3156,1.3155,1.3154,1.3153,&
-1.3152,1.3151,1.3150,1.3150,1.3149,1.3148,1.3147,1.3146,1.3146,1.3145,&
-1.3144,1.3143,1.3142,1.3141,1.3140,1.3139,1.3138,1.3137,1.3136,1.3135,&
-1.3134,1.3133,1.3132,1.3131,1.3130,1.3129,1.3128,1.3127,1.3126,1.3125,&
-1.3124,1.3123,1.3122,1.3121,1.3120,1.3119,1.3118,1.3117,1.3116,1.3115,&
-1.3114,1.3113,1.3112,1.3110,1.3109,1.3108,1.3107,1.3106,1.3105,1.3104,&
-1.3103,1.3102,1.3101,1.3099,1.3098,1.3097,1.3096,1.3095,1.3094,1.3093,&
-1.3092,1.3091,1.3090,1.3088,1.3087,1.3086,1.3085,1.3084,1.3082,1.3081,&
-1.3080,1.3079,1.3078,1.3076,1.3075,1.3074,1.3073,1.3072,1.3070,1.3069,&
-1.3068,1.3067,1.3066,1.3064,1.3063,1.3062,1.3061,1.3060,1.3059,1.3058/
- data (refractive(i),i=1701,1800)/&
-1.3057,1.3056,1.3055,1.3053,1.3052,1.3051,1.3050,1.3049,1.3047,1.3046,&
-1.3045,1.3044,1.3043,1.3041,1.3040,1.3039,1.3038,1.3037,1.3035,1.3034,&
-1.3033,1.3032,1.3031,1.3029,1.3028,1.3027,1.3026,1.3025,1.3023,1.3022,&
-1.3021,1.3020,1.3019,1.3017,1.3016,1.3015,1.3014,1.3013,1.3011,1.3010,&
-1.3009,1.3008,1.3007,1.3005,1.3004,1.3003,1.3002,1.3001,1.2999,1.2998,&
-1.2997,1.2996,1.2995,1.2993,1.2992,1.2991,1.2990,1.2989,1.2987,1.2986,&
-1.2985,1.2984,1.2983,1.2981,1.2980,1.2979,1.2978,1.2977,1.2975,1.2974,&
-1.2973,1.2972,1.2970,1.2969,1.2967,1.2966,1.2965,1.2964,1.2962,1.2961,&
-1.2960,1.2959,1.2958,1.2957,1.2956,1.2955,1.2954,1.2953,1.2951,1.2950,&
-1.2949,1.2948,1.2947,1.2945,1.2944,1.2943,1.2942,1.2941,1.2939,1.2938/
- data (refractive(i),i=1801,1900)/&
-1.2937,1.2936,1.2935,1.2933,1.2932,1.2931,1.2930,1.2928,1.2927,1.2925,&
-1.2924,1.2923,1.2922,1.2920,1.2919,1.2918,1.2917,1.2916,1.2914,1.2913,&
-1.2912,1.2911,1.2910,1.2908,1.2907,1.2906,1.2905,1.2904,1.2902,1.2901,&
-1.2900,1.2899,1.2898,1.2896,1.2895,1.2894,1.2893,1.2892,1.2890,1.2889,&
-1.2888,1.2887,1.2886,1.2884,1.2883,1.2882,1.2881,1.2880,1.2878,1.2877,&
-1.2876,1.2875,1.2874,1.2873,1.2872,1.2871,1.2870,1.2869,1.2867,1.2866,&
-1.2865,1.2864,1.2863,1.2862,1.2861,1.2860,1.2859,1.2858,1.2856,1.2855,&
-1.2854,1.2853,1.2852,1.2850,1.2849,1.2848,1.2847,1.2846,1.2845,1.2844,&
-1.2843,1.2842,1.2841,1.2840,1.2839,1.2838,1.2837,1.2836,1.2834,1.2833,&
-1.2832,1.2831,1.2830,1.2829,1.2828,1.2827,1.2826,1.2825,1.2824,1.2823/
- data (refractive(i),i=1901,2000)/&
-1.2822,1.2821,1.2820,1.2818,1.2817,1.2816,1.2815,1.2814,1.2813,1.2812,&
-1.2811,1.2810,1.2809,1.2808,1.2807,1.2806,1.2805,1.2804,1.2803,1.2802,&
-1.2801,1.2800,1.2799,1.2799,1.2798,1.2797,1.2796,1.2795,1.2795,1.2794,&
-1.2793,1.2792,1.2791,1.2791,1.2790,1.2789,1.2788,1.2788,1.2787,1.2787,&
-1.2786,1.2785,1.2785,1.2784,1.2784,1.2783,1.2782,1.2781,1.2780,1.2779,&
-1.2778,1.2777,1.2776,1.2776,1.2775,1.2774,1.2773,1.2772,1.2771,1.2770,&
-1.2769,1.2768,1.2767,1.2766,1.2765,1.2764,1.2763,1.2762,1.2761,1.2760,&
-1.2759,1.2758,1.2757,1.2757,1.2756,1.2755,1.2754,1.2753,1.2753,1.2752,&
-1.2751,1.2750,1.2749,1.2749,1.2748,1.2747,1.2746,1.2745,1.2745,1.2744,&
-1.2743,1.2742,1.2742,1.2741,1.2741,1.2740,1.2739,1.2739,1.2738,1.2738/
- data (refractive(i),i=2001,2101)/&
-1.2737,1.2736,1.2736,1.2735,1.2735,1.2734,1.2733,1.2733,1.2732,1.2732,&
-1.2731,1.2730,1.2730,1.2729,1.2729,1.2728,1.2727,1.2727,1.2726,1.2726,&
-1.2725,1.2725,1.2724,1.2724,1.2723,1.2723,1.2722,1.2722,1.2721,1.2721,&
-1.2720,1.2720,1.2719,1.2719,1.2718,1.2718,1.2717,1.2717,1.2716,1.2716,&
-1.2715,1.2715,1.2714,1.2714,1.2713,1.2713,1.2713,1.2713,1.2712,1.2712,&
-1.2712,1.2712,1.2711,1.2711,1.2710,1.2710,1.2710,1.2710,1.2709,1.2709,&
-1.2709,1.2709,1.2709,1.2708,1.2708,1.2708,1.2708,1.2708,1.2708,1.2708,&
-1.2708,1.2708,1.2708,1.2708,1.2708,1.2708,1.2708,1.2709,1.2709,1.2710,&
-1.2710,1.2711,1.2712,1.2712,1.2713,1.2714,1.2715,1.2716,1.2717,1.2718,&
-1.2719,1.2720,1.2722,1.2723,1.2725,1.2726,1.2728,1.2730,1.2732,1.2734,&
-1.2736/
-
-! ********************************************************************************
-! Specific absorption coefficient of chlorophyll
-! ********************************************************************************
-
- data (k_Cab(i),i=1,100)/&
-2.676E-02,3.113E-02,3.561E-02,3.972E-02,4.321E-02,4.606E-02,4.815E-02,5.007E-02,5.189E-02,5.346E-02,&
-5.526E-02,5.721E-02,5.891E-02,6.004E-02,6.100E-02,6.189E-02,6.259E-02,6.318E-02,6.367E-02,6.416E-02,&
-6.475E-02,6.555E-02,6.628E-02,6.690E-02,6.718E-02,6.747E-02,6.795E-02,6.855E-02,6.897E-02,6.931E-02,&
-6.979E-02,7.037E-02,7.082E-02,7.114E-02,7.141E-02,7.160E-02,7.174E-02,7.177E-02,7.171E-02,7.156E-02,&
-7.120E-02,7.048E-02,6.963E-02,6.883E-02,6.789E-02,6.685E-02,6.583E-02,6.480E-02,6.358E-02,6.225E-02,&
-6.090E-02,5.956E-02,5.813E-02,5.677E-02,5.553E-02,5.439E-02,5.330E-02,5.224E-02,5.135E-02,5.056E-02,&
-4.990E-02,4.934E-02,4.891E-02,4.852E-02,4.812E-02,4.767E-02,4.729E-02,4.698E-02,4.673E-02,4.647E-02,&
-4.622E-02,4.606E-02,4.589E-02,4.570E-02,4.552E-02,4.535E-02,4.519E-02,4.492E-02,4.466E-02,4.442E-02,&
-4.413E-02,4.373E-02,4.325E-02,4.275E-02,4.216E-02,4.149E-02,4.073E-02,3.989E-02,3.895E-02,3.791E-02,&
-3.683E-02,3.571E-02,3.455E-02,3.335E-02,3.215E-02,3.094E-02,2.972E-02,2.850E-02,2.730E-02,2.615E-02/
- data (k_Cab(i),i=101,200)/&
-2.502e-02,2.392e-02,2.284e-02,2.179e-02,2.077e-02,1.979e-02,1.883e-02,1.790e-02,1.702e-02,1.617e-02,&
-1.536e-02,1.458e-02,1.383e-02,1.312e-02,1.246e-02,1.183e-02,1.126e-02,1.074e-02,1.027e-02,9.855e-03,&
-9.490e-03,9.176e-03,8.910e-03,8.698e-03,8.542e-03,8.437e-03,8.381e-03,8.373e-03,8.411e-03,8.486e-03,&
-8.595e-03,8.735e-03,8.905e-03,9.098e-03,9.313e-03,9.544e-03,9.789e-03,1.004e-02,1.030e-02,1.056e-02,&
-1.083e-02,1.108e-02,1.133e-02,1.157e-02,1.180e-02,1.201e-02,1.221e-02,1.239e-02,1.256e-02,1.271e-02,&
-1.286e-02,1.300e-02,1.314e-02,1.327e-02,1.340e-02,1.353e-02,1.367e-02,1.380e-02,1.394e-02,1.408e-02,&
-1.421e-02,1.435e-02,1.451e-02,1.469e-02,1.489e-02,1.511e-02,1.535e-02,1.561e-02,1.589e-02,1.619e-02,&
-1.650e-02,1.683e-02,1.716e-02,1.749e-02,1.782e-02,1.814e-02,1.846e-02,1.876e-02,1.905e-02,1.933e-02,&
-1.960e-02,1.986e-02,2.011e-02,2.035e-02,2.058e-02,2.079e-02,2.100e-02,2.120e-02,2.138e-02,2.155e-02,&
-2.171e-02,2.185e-02,2.199e-02,2.211e-02,2.222e-02,2.231e-02,2.241e-02,2.250e-02,2.259e-02,2.269e-02/
- data (k_Cab(i),i=201,300)/&
-2.279e-02,2.291e-02,2.304e-02,2.319e-02,2.337e-02,2.356e-02,2.378e-02,2.403e-02,2.430e-02,2.459e-02,&
-2.489e-02,2.520e-02,2.552e-02,2.586e-02,2.620e-02,2.654e-02,2.687e-02,2.719e-02,2.750e-02,2.778e-02,&
-2.803e-02,2.826e-02,2.845e-02,2.862e-02,2.877e-02,2.888e-02,2.898e-02,2.905e-02,2.910e-02,2.916e-02,&
-2.921e-02,2.928e-02,2.938e-02,2.951e-02,2.969e-02,2.993e-02,3.022e-02,3.059e-02,3.102e-02,3.152e-02,&
-3.208e-02,3.269e-02,3.335e-02,3.405e-02,3.478e-02,3.552e-02,3.626e-02,3.698e-02,3.767e-02,3.831e-02,&
-3.890e-02,3.943e-02,3.993e-02,4.043e-02,4.095e-02,4.152e-02,4.217e-02,4.293e-02,4.380e-02,4.478e-02,&
-4.588e-02,4.709e-02,4.837e-02,4.970e-02,5.106e-02,5.239e-02,5.368e-02,5.488e-02,5.595e-02,5.688e-02,&
-5.765e-02,5.826e-02,5.871e-02,5.900e-02,5.914e-02,5.914e-02,5.899e-02,5.868e-02,5.819e-02,5.748e-02,&
-5.650e-02,5.522e-02,5.361e-02,5.168e-02,4.941e-02,4.689e-02,4.416e-02,4.128e-02,3.834e-02,3.541e-02,&
-3.254e-02,2.978e-02,2.717e-02,2.473e-02,2.247e-02,2.039e-02,1.850e-02,1.678e-02,1.523e-02,1.384e-02/
- data (k_Cab(i),i=301,400)/&
-1.259e-02,1.147e-02,1.048e-02,9.594e-03,8.804e-03,8.099e-03,7.465e-03,6.892e-03,6.371e-03,5.897e-03,&
-5.461e-03,5.061e-03,4.691e-03,4.348e-03,4.029e-03,3.733e-03,3.456e-03,3.199e-03,2.960e-03,2.736e-03,&
-2.528e-03,2.334e-03,2.153e-03,1.985e-03,1.829e-03,1.683e-03,1.549e-03,1.424e-03,1.309e-03,1.202e-03,&
-1.103e-03,1.012e-03,9.269e-04,8.485e-04,7.768e-04,7.110e-04,6.505e-04,5.948e-04,5.437e-04,4.967e-04,&
-4.539e-04,4.150e-04,3.794e-04,3.471e-04,3.177e-04,2.909e-04,2.666e-04,2.450e-04,2.255e-04,2.079e-04,&
-1.924e-04,1.773e-04,1.629e-04,1.493e-04,1.364e-04,1.243e-04,1.128e-04,1.020e-04,9.188e-05,8.240e-05,&
-7.354e-05,6.529e-05,5.763e-05,5.054e-05,4.400e-05,3.800e-05,3.251e-05,2.753e-05,2.302e-05,1.898e-05,&
-1.539e-05,1.222e-05,9.466e-06,7.101e-06,5.109e-06,3.473e-06,2.175e-06,1.196e-06,5.200e-07,1.270e-07,&
-0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,&
-0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
- data (k_Cab(i),i=401,2101)/1701*0./
-
-! ********************************************************************************
-! Specific absorption coefficient of carotenoids
-! ********************************************************************************
-
- data (k_Car(i),i=1,100)/&
-2.895e-01,2.796e-01,2.693e-01,2.601e-01,2.525e-01,2.456e-01,2.402e-01,2.352e-01,2.303e-01,2.255e-01,&
-2.203e-01,2.151e-01,2.105e-01,2.069e-01,2.035e-01,2.002e-01,1.971e-01,1.943e-01,1.914e-01,1.889e-01,&
-1.866e-01,1.844e-01,1.823e-01,1.805e-01,1.792e-01,1.777e-01,1.762e-01,1.747e-01,1.735e-01,1.724e-01,&
-1.715e-01,1.707e-01,1.701e-01,1.694e-01,1.688e-01,1.682e-01,1.678e-01,1.674e-01,1.670e-01,1.667e-01,&
-1.664e-01,1.659e-01,1.653e-01,1.647e-01,1.641e-01,1.632e-01,1.624e-01,1.615e-01,1.605e-01,1.592e-01,&
-1.577e-01,1.561e-01,1.544e-01,1.528e-01,1.513e-01,1.497e-01,1.483e-01,1.468e-01,1.454e-01,1.440e-01,&
-1.425e-01,1.411e-01,1.397e-01,1.384e-01,1.371e-01,1.358e-01,1.347e-01,1.335e-01,1.325e-01,1.316e-01,&
-1.308e-01,1.302e-01,1.296e-01,1.292e-01,1.288e-01,1.285e-01,1.282e-01,1.280e-01,1.277e-01,1.274e-01,&
-1.270e-01,1.265e-01,1.260e-01,1.253e-01,1.244e-01,1.234e-01,1.223e-01,1.211e-01,1.197e-01,1.182e-01,&
-1.166e-01,1.150e-01,1.132e-01,1.114e-01,1.096e-01,1.076e-01,1.056e-01,1.035e-01,1.013e-01,9.901e-02/
- data (k_Car(i),i=101,200)/&
-9.660e-02,9.420e-02,9.179e-02,8.939e-02,8.699e-02,8.460e-02,8.221e-02,7.983e-02,7.746e-02,7.510e-02,&
-7.275e-02,7.041e-02,6.809e-02,6.578e-02,6.349e-02,6.122e-02,5.897e-02,5.673e-02,5.453e-02,5.234e-02,&
-5.018e-02,4.805e-02,4.594e-02,4.387e-02,4.182e-02,3.981e-02,3.783e-02,3.588e-02,3.397e-02,3.210e-02,&
-3.026e-02,2.847e-02,2.672e-02,2.501e-02,2.334e-02,2.172e-02,2.015e-02,1.863e-02,1.715e-02,1.573e-02,&
-1.436e-02,1.304e-02,1.178e-02,1.057e-02,9.420e-03,8.330e-03,7.301e-03,6.334e-03,5.429e-03,4.589e-03,&
-3.816e-03,3.109e-03,2.471e-03,1.903e-03,1.406e-03,9.822e-04,6.322e-04,3.577e-04,1.599e-04,4.019e-05,&
-0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,&
-0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,&
-0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,&
-0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00,0.000e+00/
- data (k_Car(i),i=201,2101)/1901*0./
-
-! ********************************************************************************
-! Specific absorption coefficient of brown pigments
-! ********************************************************************************
-
- data (k_Brown(i),i=1,100)/&
-5.272e-01,5.262e-01,5.252e-01,5.242e-01,5.232e-01,5.222e-01,5.212e-01,5.202e-01,5.192e-01,5.182e-01,&
-5.172e-01,5.162e-01,5.152e-01,5.142e-01,5.132e-01,5.122e-01,5.112e-01,5.102e-01,5.092e-01,5.082e-01,&
-5.072e-01,5.062e-01,5.052e-01,5.042e-01,5.032e-01,5.022e-01,5.012e-01,5.002e-01,4.992e-01,4.982e-01,&
-4.972e-01,4.960e-01,4.948e-01,4.936e-01,4.924e-01,4.912e-01,4.900e-01,4.888e-01,4.876e-01,4.864e-01,&
-4.852e-01,4.840e-01,4.829e-01,4.817e-01,4.805e-01,4.793e-01,4.781e-01,4.769e-01,4.757e-01,4.745e-01,&
-4.733e-01,4.720e-01,4.708e-01,4.695e-01,4.683e-01,4.670e-01,4.658e-01,4.645e-01,4.633e-01,4.620e-01,&
-4.608e-01,4.597e-01,4.587e-01,4.577e-01,4.566e-01,4.556e-01,4.546e-01,4.535e-01,4.525e-01,4.515e-01,&
-4.504e-01,4.494e-01,4.484e-01,4.473e-01,4.463e-01,4.453e-01,4.442e-01,4.432e-01,4.422e-01,4.411e-01,&
-4.401e-01,4.388e-01,4.375e-01,4.362e-01,4.350e-01,4.337e-01,4.324e-01,4.311e-01,4.298e-01,4.285e-01,&
-4.272e-01,4.260e-01,4.247e-01,4.234e-01,4.221e-01,4.208e-01,4.195e-01,4.183e-01,4.170e-01,4.157e-01/
- data (k_Brown(i),i=101,200)/&
-4.144e-01,4.135e-01,4.127e-01,4.118e-01,4.109e-01,4.100e-01,4.092e-01,4.083e-01,4.074e-01,4.065e-01,&
-4.057e-01,4.047e-01,4.038e-01,4.028e-01,4.019e-01,4.009e-01,4.000e-01,3.990e-01,3.981e-01,3.971e-01,&
-3.962e-01,3.952e-01,3.943e-01,3.933e-01,3.924e-01,3.914e-01,3.905e-01,3.895e-01,3.886e-01,3.876e-01,&
-3.867e-01,3.846e-01,3.824e-01,3.803e-01,3.782e-01,3.760e-01,3.739e-01,3.718e-01,3.696e-01,3.675e-01,&
-3.654e-01,3.639e-01,3.625e-01,3.611e-01,3.597e-01,3.582e-01,3.568e-01,3.554e-01,3.540e-01,3.525e-01,&
-3.511e-01,3.500e-01,3.489e-01,3.478e-01,3.467e-01,3.456e-01,3.445e-01,3.434e-01,3.423e-01,3.412e-01,&
-3.401e-01,3.383e-01,3.366e-01,3.348e-01,3.330e-01,3.312e-01,3.294e-01,3.276e-01,3.258e-01,3.241e-01,&
-3.223e-01,3.202e-01,3.182e-01,3.162e-01,3.141e-01,3.121e-01,3.100e-01,3.080e-01,3.059e-01,3.039e-01,&
-3.019e-01,2.999e-01,2.979e-01,2.959e-01,2.940e-01,2.920e-01,2.900e-01,2.881e-01,2.861e-01,2.841e-01,&
-2.821e-01,2.803e-01,2.784e-01,2.766e-01,2.747e-01,2.728e-01,2.710e-01,2.691e-01,2.673e-01,2.654e-01/
- data (k_Brown(i),i=201,300)/&
-2.636e-01,2.618e-01,2.601e-01,2.584e-01,2.566e-01,2.549e-01,2.532e-01,2.515e-01,2.497e-01,2.480e-01,&
-2.463e-01,2.447e-01,2.431e-01,2.414e-01,2.398e-01,2.382e-01,2.366e-01,2.350e-01,2.334e-01,2.318e-01,&
-2.302e-01,2.288e-01,2.273e-01,2.258e-01,2.244e-01,2.229e-01,2.215e-01,2.200e-01,2.185e-01,2.171e-01,&
-2.156e-01,2.143e-01,2.129e-01,2.115e-01,2.102e-01,2.088e-01,2.074e-01,2.061e-01,2.047e-01,2.033e-01,&
-2.020e-01,2.007e-01,1.994e-01,1.981e-01,1.968e-01,1.955e-01,1.942e-01,1.929e-01,1.916e-01,1.903e-01,&
-1.890e-01,1.878e-01,1.865e-01,1.853e-01,1.841e-01,1.829e-01,1.816e-01,1.804e-01,1.792e-01,1.780e-01,&
-1.768e-01,1.755e-01,1.742e-01,1.729e-01,1.717e-01,1.704e-01,1.691e-01,1.679e-01,1.666e-01,1.653e-01,&
-1.641e-01,1.627e-01,1.613e-01,1.600e-01,1.586e-01,1.572e-01,1.559e-01,1.545e-01,1.532e-01,1.518e-01,&
-1.504e-01,1.491e-01,1.478e-01,1.464e-01,1.451e-01,1.437e-01,1.424e-01,1.411e-01,1.397e-01,1.384e-01,&
-1.370e-01,1.358e-01,1.345e-01,1.333e-01,1.320e-01,1.308e-01,1.295e-01,1.283e-01,1.270e-01,1.258e-01/
- data (k_Brown(i),i=301,400)/&
-1.245e-01,1.234e-01,1.223e-01,1.212e-01,1.200e-01,1.189e-01,1.178e-01,1.167e-01,1.156e-01,1.144e-01,&
-1.133e-01,1.122e-01,1.111e-01,1.100e-01,1.089e-01,1.078e-01,1.067e-01,1.056e-01,1.046e-01,1.035e-01,&
-1.024e-01,1.013e-01,1.003e-01,9.931e-02,9.829e-02,9.727e-02,9.625e-02,9.524e-02,9.422e-02,9.320e-02,&
-9.218e-02,9.120e-02,9.022e-02,8.924e-02,8.827e-02,8.729e-02,8.631e-02,8.533e-02,8.435e-02,8.337e-02,&
-8.239e-02,8.155e-02,8.070e-02,7.985e-02,7.901e-02,7.816e-02,7.732e-02,7.647e-02,7.562e-02,7.478e-02,&
-7.393e-02,7.319e-02,7.245e-02,7.171e-02,7.097e-02,7.023e-02,6.949e-02,6.875e-02,6.801e-02,6.727e-02,&
-6.653e-02,6.586e-02,6.519e-02,6.452e-02,6.385e-02,6.318e-02,6.251e-02,6.184e-02,6.117e-02,6.050e-02,&
-5.983e-02,5.913e-02,5.843e-02,5.773e-02,5.704e-02,5.634e-02,5.564e-02,5.494e-02,5.424e-02,5.354e-02,&
-5.284e-02,5.226e-02,5.167e-02,5.109e-02,5.050e-02,4.992e-02,4.933e-02,4.874e-02,4.816e-02,4.757e-02,&
-4.699e-02,4.646e-02,4.594e-02,4.542e-02,4.490e-02,4.437e-02,4.385e-02,4.333e-02,4.281e-02,4.228e-02/
- data (k_Brown(i),i=401,500)/&
-4.176e-02,4.128e-02,4.081e-02,4.033e-02,3.985e-02,3.937e-02,3.889e-02,3.841e-02,3.793e-02,3.746e-02,&
-3.698e-02,3.657e-02,3.615e-02,3.574e-02,3.533e-02,3.492e-02,3.451e-02,3.409e-02,3.368e-02,3.327e-02,&
-3.286e-02,3.249e-02,3.212e-02,3.175e-02,3.139e-02,3.102e-02,3.065e-02,3.028e-02,2.991e-02,2.954e-02,&
-2.918e-02,2.885e-02,2.852e-02,2.819e-02,2.786e-02,2.753e-02,2.720e-02,2.688e-02,2.655e-02,2.622e-02,&
-2.589e-02,2.559e-02,2.529e-02,2.499e-02,2.469e-02,2.440e-02,2.410e-02,2.380e-02,2.350e-02,2.320e-02,&
-2.290e-02,2.264e-02,2.238e-02,2.212e-02,2.186e-02,2.159e-02,2.133e-02,2.107e-02,2.081e-02,2.055e-02,&
-2.029e-02,2.006e-02,1.983e-02,1.961e-02,1.938e-02,1.915e-02,1.893e-02,1.870e-02,1.847e-02,1.825e-02,&
-1.802e-02,1.782e-02,1.762e-02,1.742e-02,1.723e-02,1.703e-02,1.683e-02,1.663e-02,1.643e-02,1.623e-02,&
-1.604e-02,1.586e-02,1.568e-02,1.551e-02,1.533e-02,1.516e-02,1.498e-02,1.481e-02,1.463e-02,1.446e-02,&
-1.428e-02,1.414e-02,1.400e-02,1.385e-02,1.371e-02,1.357e-02,1.343e-02,1.328e-02,1.314e-02,1.300e-02/
- data (k_Brown(i),i=501,600)/&
-1.286e-02,1.275e-02,1.265e-02,1.255e-02,1.245e-02,1.235e-02,1.225e-02,1.215e-02,1.205e-02,1.195e-02,&
-1.185e-02,1.175e-02,1.165e-02,1.155e-02,1.145e-02,1.135e-02,1.125e-02,1.115e-02,1.105e-02,1.096e-02,&
-1.086e-02,1.076e-02,1.066e-02,1.056e-02,1.047e-02,1.037e-02,1.027e-02,1.017e-02,1.008e-02,9.980e-03,&
-9.884e-03,9.787e-03,9.691e-03,9.595e-03,9.500e-03,9.404e-03,9.309e-03,9.214e-03,9.120e-03,9.025e-03,&
-8.931e-03,8.837e-03,8.743e-03,8.650e-03,8.557e-03,8.464e-03,8.371e-03,8.279e-03,8.187e-03,8.095e-03,&
-8.004e-03,7.913e-03,7.822e-03,7.732e-03,7.641e-03,7.552e-03,7.462e-03,7.373e-03,7.284e-03,7.195e-03,&
-7.107e-03,7.019e-03,6.932e-03,6.844e-03,6.758e-03,6.671e-03,6.585e-03,6.499e-03,6.414e-03,6.329e-03,&
-6.244e-03,6.160e-03,6.076e-03,5.993e-03,5.910e-03,5.827e-03,5.745e-03,5.663e-03,5.581e-03,5.500e-03,&
-5.419e-03,5.339e-03,5.259e-03,5.180e-03,5.101e-03,5.023e-03,4.945e-03,4.867e-03,4.790e-03,4.713e-03,&
-4.637e-03,4.561e-03,4.486e-03,4.411e-03,4.337e-03,4.263e-03,4.190e-03,4.117e-03,4.044e-03,3.972e-03/
- data (k_Brown(i),i=601,700)/&
-3.901e-03,3.830e-03,3.760e-03,3.690e-03,3.620e-03,3.552e-03,3.483e-03,3.416e-03,3.348e-03,3.282e-03,&
-3.215e-03,3.150e-03,3.085e-03,3.020e-03,2.956e-03,2.893e-03,2.830e-03,2.768e-03,2.706e-03,2.645e-03,&
-2.585e-03,2.525e-03,2.465e-03,2.407e-03,2.348e-03,2.291e-03,2.234e-03,2.178e-03,2.122e-03,2.067e-03,&
-2.012e-03,1.959e-03,1.905e-03,1.853e-03,1.801e-03,1.750e-03,1.699e-03,1.649e-03,1.600e-03,1.551e-03,&
-1.503e-03,1.456e-03,1.409e-03,1.363e-03,1.318e-03,1.274e-03,1.230e-03,1.187e-03,1.144e-03,1.102e-03,&
-1.061e-03,1.021e-03,9.811e-04,9.422e-04,9.040e-04,8.665e-04,8.297e-04,7.937e-04,7.584e-04,7.239e-04,&
-6.901e-04,6.571e-04,6.248e-04,5.933e-04,5.626e-04,5.326e-04,5.034e-04,4.750e-04,4.473e-04,4.205e-04,&
-3.944e-04,3.691e-04,3.446e-04,3.210e-04,2.981e-04,2.760e-04,2.548e-04,2.344e-04,2.148e-04,1.960e-04,&
-1.780e-04,1.609e-04,1.446e-04,1.292e-04,1.146e-04,1.009e-04,8.804e-05,7.603e-05,6.488e-05,5.460e-05,&
-4.519e-05,3.666e-05,2.901e-05,2.225e-05,1.637e-05,1.138e-05,7.297e-06,4.111e-06,1.830e-06,4.581e-07/
- data (k_Brown(i),i=701,2101)/1401*0./
-
-! ********************************************************************************
-! Specific absorption coefficient of water
-! ********************************************************************************
-
- data (k_Cw(i),i=1,100)/&
-5.800E-05,5.852E-05,5.900E-05,5.989E-05,6.100E-05,6.203E-05,6.300E-05,6.399E-05,6.500E-05,6.603E-05,&
-6.700E-05,6.790E-05,6.900E-05,7.050E-05,7.200E-05,7.312E-05,7.400E-05,7.490E-05,7.600E-05,7.740E-05,&
-7.900E-05,8.063E-05,8.200E-05,8.297E-05,8.400E-05,8.551E-05,8.700E-05,8.800E-05,8.900E-05,9.050E-05,&
-9.200E-05,9.300E-05,9.400E-05,9.550E-05,9.700E-05,9.801E-05,9.900E-05,1.005E-04,1.020E-04,1.031E-04,&
-1.040E-04,1.050E-04,1.060E-04,1.070E-04,1.080E-04,1.090E-04,1.100E-04,1.110E-04,1.120E-04,1.130E-04,&
-1.140E-04,1.150E-04,1.160E-04,1.170E-04,1.180E-04,1.190E-04,1.200E-04,1.210E-04,1.220E-04,1.230E-04,&
-1.240E-04,1.250E-04,1.260E-04,1.270E-04,1.280E-04,1.289E-04,1.300E-04,1.315E-04,1.330E-04,1.340E-04,&
-1.350E-04,1.364E-04,1.380E-04,1.396E-04,1.410E-04,1.424E-04,1.440E-04,1.459E-04,1.480E-04,1.499E-04,&
-1.520E-04,1.544E-04,1.570E-04,1.596E-04,1.620E-04,1.643E-04,1.670E-04,1.704E-04,1.740E-04,1.775E-04,&
-1.810E-04,1.849E-04,1.890E-04,1.934E-04,1.980E-04,2.031E-04,2.090E-04,2.158E-04,2.230E-04,2.303E-04/
- data (k_Cw(i),i=101,200)/&
-2.380E-04,2.463E-04,2.550E-04,2.640E-04,2.730E-04,2.819E-04,2.910E-04,3.004E-04,3.100E-04,3.194E-04,&
-3.290E-04,3.390E-04,3.490E-04,3.588E-04,3.680E-04,3.767E-04,3.860E-04,3.962E-04,4.040E-04,4.069E-04,&
-4.090E-04,4.138E-04,4.160E-04,4.112E-04,4.090E-04,4.176E-04,4.270E-04,4.268E-04,4.230E-04,4.237E-04,&
-4.290E-04,4.371E-04,4.450E-04,4.506E-04,4.560E-04,4.631E-04,4.700E-04,4.748E-04,4.800E-04,4.879E-04,&
-4.950E-04,4.983E-04,5.030E-04,5.141E-04,5.270E-04,5.363E-04,5.440E-04,5.532E-04,5.640E-04,5.759E-04,&
-5.880E-04,5.998E-04,6.110E-04,6.215E-04,6.310E-04,6.391E-04,6.460E-04,6.520E-04,6.580E-04,6.647E-04,&
-6.720E-04,6.793E-04,6.860E-04,6.920E-04,6.990E-04,7.084E-04,7.180E-04,7.257E-04,7.340E-04,7.455E-04,&
-7.590E-04,7.729E-04,7.870E-04,8.020E-04,8.190E-04,8.386E-04,8.580E-04,8.754E-04,8.960E-04,9.238E-04,&
-9.520E-04,9.745E-04,1.000E-03,1.037E-03,1.079E-03,1.119E-03,1.159E-03,1.204E-03,1.253E-03,1.304E-03,&
-1.356E-03,1.408E-03,1.459E-03,1.510E-03,1.567E-03,1.635E-03,1.700E-03,1.758E-03,1.860E-03,2.042E-03/
- data (k_Cw(i),i=201,300)/&
-2.224E-03,2.323E-03,2.366E-03,2.400E-03,2.448E-03,2.519E-03,2.587E-03,2.629E-03,2.653E-03,2.674E-03,&
-2.691E-03,2.704E-03,2.715E-03,2.727E-03,2.740E-03,2.753E-03,2.764E-03,2.775E-03,2.785E-03,2.797E-03,&
-2.810E-03,2.824E-03,2.839E-03,2.854E-03,2.868E-03,2.881E-03,2.893E-03,2.907E-03,2.922E-03,2.938E-03,&
-2.955E-03,2.972E-03,2.988E-03,3.000E-03,3.011E-03,3.023E-03,3.038E-03,3.057E-03,3.076E-03,3.094E-03,&
-3.111E-03,3.127E-03,3.144E-03,3.162E-03,3.181E-03,3.202E-03,3.223E-03,3.242E-03,3.263E-03,3.289E-03,&
-3.315E-03,3.338E-03,3.362E-03,3.390E-03,3.423E-03,3.461E-03,3.508E-03,3.567E-03,3.636E-03,3.712E-03,&
-3.791E-03,3.866E-03,3.931E-03,3.981E-03,4.019E-03,4.049E-03,4.072E-03,4.087E-03,4.098E-03,4.109E-03,&
-4.122E-03,4.137E-03,4.150E-03,4.160E-03,4.173E-03,4.196E-03,4.223E-03,4.248E-03,4.270E-03,4.293E-03,&
-4.318E-03,4.347E-03,4.381E-03,4.418E-03,4.458E-03,4.500E-03,4.545E-03,4.594E-03,4.646E-03,4.701E-03,&
-4.760E-03,4.827E-03,4.903E-03,4.986E-03,5.071E-03,5.154E-03,5.244E-03,5.351E-03,5.470E-03,5.594E-03/
- data (k_Cw(i),i=301,400)/&
-5.722E-03,5.855E-03,5.995E-03,6.146E-03,6.303E-03,6.463E-03,6.628E-03,6.804E-03,6.993E-03,7.197E-03,&
-7.415E-03,7.647E-03,7.893E-03,8.157E-03,8.445E-03,8.763E-03,9.109E-03,9.479E-03,9.871E-03,1.029E-02,&
-1.072E-02,1.119E-02,1.168E-02,1.218E-02,1.268E-02,1.319E-02,1.372E-02,1.428E-02,1.487E-02,1.551E-02,&
-1.621E-02,1.699E-02,1.787E-02,1.886E-02,1.992E-02,2.101E-02,2.207E-02,2.306E-02,2.394E-02,2.469E-02,&
-2.532E-02,2.583E-02,2.623E-02,2.652E-02,2.672E-02,2.689E-02,2.702E-02,2.713E-02,2.722E-02,2.728E-02,&
-2.733E-02,2.738E-02,2.741E-02,2.745E-02,2.748E-02,2.751E-02,2.754E-02,2.758E-02,2.763E-02,2.767E-02,&
-2.771E-02,2.773E-02,2.773E-02,2.774E-02,2.774E-02,2.773E-02,2.770E-02,2.766E-02,2.761E-02,2.757E-02,&
-2.754E-02,2.752E-02,2.748E-02,2.741E-02,2.731E-02,2.720E-02,2.710E-02,2.701E-02,2.690E-02,2.675E-02,&
-2.659E-02,2.645E-02,2.633E-02,2.624E-02,2.613E-02,2.593E-02,2.558E-02,2.523E-02,2.513E-02,2.501E-02,&
-2.466E-02,2.447E-02,2.412E-02,2.389E-02,2.374E-02,2.355E-02,2.337E-02,2.318E-02,2.304E-02,2.281E-02/
- data (k_Cw(i),i=401,500)/&
-2.246E-02,2.243E-02,2.238E-02,2.222E-02,2.204E-02,2.201E-02,2.204E-02,2.196E-02,2.177E-02,2.190E-02,&
-2.188E-02,2.188E-02,2.198E-02,2.210E-02,2.223E-02,2.233E-02,2.248E-02,2.276E-02,2.304E-02,2.311E-02,&
-2.329E-02,2.388E-02,2.446E-02,2.475E-02,2.516E-02,2.620E-02,2.769E-02,2.830E-02,2.914E-02,3.108E-02,&
-3.214E-02,3.297E-02,3.459E-02,3.606E-02,3.662E-02,3.702E-02,3.788E-02,3.829E-02,3.854E-02,3.909E-02,&
-3.949E-02,3.972E-02,4.000E-02,4.040E-02,4.057E-02,4.075E-02,4.115E-02,4.127E-02,4.149E-02,4.204E-02,&
-4.199E-02,4.223E-02,4.254E-02,4.272E-02,4.280E-02,4.306E-02,4.360E-02,4.369E-02,4.379E-02,4.433E-02,&
-4.454E-02,4.466E-02,4.505E-02,4.527E-02,4.552E-02,4.605E-02,4.658E-02,4.691E-02,4.705E-02,4.713E-02,&
-4.752E-02,4.833E-02,4.867E-02,4.894E-02,4.960E-02,5.006E-02,5.050E-02,5.115E-02,5.153E-02,5.204E-02,&
-5.298E-02,5.346E-02,5.386E-02,5.465E-02,5.528E-02,5.566E-02,5.596E-02,5.653E-02,5.745E-02,5.789E-02,&
-5.831E-02,5.924E-02,5.982E-02,6.009E-02,6.035E-02,6.094E-02,6.185E-02,6.226E-02,6.269E-02,6.360E-02/
- data (k_Cw(i),i=501,600)/&
-6.407E-02,6.458E-02,6.562E-02,6.636E-02,6.672E-02,6.699E-02,6.769E-02,6.900E-02,6.989E-02,7.037E-02,&
-7.085E-02,7.187E-02,7.358E-02,7.486E-02,7.562E-02,7.630E-02,7.792E-02,8.085E-02,8.292E-02,8.410E-02,&
-8.528E-02,8.801E-02,9.268E-02,9.584E-02,9.819E-02,1.012E-01,1.042E-01,1.066E-01,1.113E-01,1.194E-01,&
-1.246E-01,1.281E-01,1.327E-01,1.374E-01,1.410E-01,1.465E-01,1.557E-01,1.635E-01,1.688E-01,1.732E-01,&
-1.818E-01,1.963E-01,2.050E-01,2.106E-01,2.187E-01,2.287E-01,2.386E-01,2.468E-01,2.542E-01,2.701E-01,&
-2.976E-01,3.153E-01,3.274E-01,3.438E-01,3.622E-01,3.785E-01,3.930E-01,4.068E-01,4.184E-01,4.273E-01,&
-4.385E-01,4.538E-01,4.611E-01,4.633E-01,4.663E-01,4.701E-01,4.733E-01,4.756E-01,4.772E-01,4.785E-01,&
-4.800E-01,4.814E-01,4.827E-01,4.843E-01,4.864E-01,4.870E-01,4.867E-01,4.864E-01,4.857E-01,4.841E-01,&
-4.821E-01,4.804E-01,4.786E-01,4.764E-01,4.738E-01,4.710E-01,4.677E-01,4.641E-01,4.604E-01,4.570E-01,&
-4.532E-01,4.482E-01,4.434E-01,4.397E-01,4.362E-01,4.316E-01,4.265E-01,4.215E-01,4.168E-01,4.121E-01/
- data (k_Cw(i),i=601,700)/&
-4.072E-01,4.017E-01,3.963E-01,3.915E-01,3.868E-01,3.816E-01,3.760E-01,3.701E-01,3.640E-01,3.581E-01,&
-3.521E-01,3.461E-01,3.402E-01,3.349E-01,3.297E-01,3.243E-01,3.191E-01,3.141E-01,3.086E-01,3.022E-01,&
-2.957E-01,2.897E-01,2.840E-01,2.782E-01,2.724E-01,2.672E-01,2.621E-01,2.566E-01,2.506E-01,2.444E-01,&
-2.391E-01,2.357E-01,2.331E-01,2.299E-01,2.251E-01,2.198E-01,2.151E-01,2.109E-01,2.064E-01,2.020E-01,&
-1.981E-01,1.944E-01,1.904E-01,1.868E-01,1.841E-01,1.818E-01,1.790E-01,1.752E-01,1.715E-01,1.687E-01,&
-1.664E-01,1.639E-01,1.613E-01,1.586E-01,1.562E-01,1.545E-01,1.532E-01,1.522E-01,1.510E-01,1.495E-01,&
-1.475E-01,1.457E-01,1.447E-01,1.442E-01,1.438E-01,1.433E-01,1.426E-01,1.418E-01,1.412E-01,1.410E-01,&
-1.409E-01,1.408E-01,1.406E-01,1.405E-01,1.408E-01,1.414E-01,1.426E-01,1.435E-01,1.438E-01,1.439E-01,&
-1.443E-01,1.456E-01,1.475E-01,1.498E-01,1.519E-01,1.534E-01,1.547E-01,1.561E-01,1.580E-01,1.604E-01,&
-1.632E-01,1.659E-01,1.677E-01,1.693E-01,1.712E-01,1.739E-01,1.777E-01,1.824E-01,1.866E-01,1.890E-01/
- data (k_Cw(i),i=701,800)/&
-1.906E-01,1.929E-01,1.967E-01,2.005E-01,2.031E-01,2.051E-01,2.079E-01,2.123E-01,2.166E-01,2.196E-01,&
-2.219E-01,2.251E-01,2.298E-01,2.337E-01,2.346E-01,2.342E-01,2.353E-01,2.397E-01,2.450E-01,2.491E-01,&
-2.528E-01,2.578E-01,2.650E-01,2.719E-01,2.765E-01,2.811E-01,2.891E-01,3.023E-01,3.164E-01,3.271E-01,&
-3.378E-01,3.533E-01,3.770E-01,4.037E-01,4.281E-01,4.502E-01,4.712E-01,4.932E-01,5.202E-01,5.572E-01,&
-6.052E-01,6.520E-01,6.863E-01,7.159E-01,7.535E-01,8.064E-01,8.597E-01,8.981E-01,9.253E-01,9.493E-01,&
-9.769E-01,1.008E+00,1.041E+00,1.073E+00,1.100E+00,1.119E+00,1.131E+00,1.140E+00,1.150E+00,1.160E+00,&
-1.170E+00,1.181E+00,1.190E+00,1.194E+00,1.196E+00,1.197E+00,1.200E+00,1.203E+00,1.205E+00,1.206E+00,&
-1.207E+00,1.213E+00,1.223E+00,1.232E+00,1.234E+00,1.232E+00,1.229E+00,1.230E+00,1.233E+00,1.236E+00,&
-1.239E+00,1.241E+00,1.244E+00,1.248E+00,1.252E+00,1.256E+00,1.258E+00,1.260E+00,1.262E+00,1.265E+00,&
-1.267E+00,1.270E+00,1.272E+00,1.275E+00,1.277E+00,1.280E+00,1.282E+00,1.283E+00,1.283E+00,1.279E+00/
- data (k_Cw(i),i=801,900)/&
-1.272E+00,1.266E+00,1.267E+00,1.271E+00,1.273E+00,1.271E+00,1.265E+00,1.260E+00,1.258E+00,1.258E+00,&
-1.257E+00,1.252E+00,1.247E+00,1.243E+00,1.243E+00,1.243E+00,1.240E+00,1.233E+00,1.224E+00,1.216E+00,&
-1.214E+00,1.214E+00,1.213E+00,1.210E+00,1.205E+00,1.200E+00,1.199E+00,1.198E+00,1.197E+00,1.194E+00,&
-1.189E+00,1.184E+00,1.180E+00,1.176E+00,1.171E+00,1.166E+00,1.161E+00,1.158E+00,1.157E+00,1.157E+00,&
-1.155E+00,1.152E+00,1.148E+00,1.142E+00,1.138E+00,1.133E+00,1.130E+00,1.126E+00,1.123E+00,1.120E+00,&
-1.116E+00,1.111E+00,1.107E+00,1.103E+00,1.101E+00,1.101E+00,1.101E+00,1.101E+00,1.100E+00,1.098E+00,&
-1.094E+00,1.089E+00,1.085E+00,1.084E+00,1.083E+00,1.083E+00,1.082E+00,1.081E+00,1.080E+00,1.079E+00,&
-1.080E+00,1.083E+00,1.087E+00,1.093E+00,1.099E+00,1.104E+00,1.107E+00,1.109E+00,1.111E+00,1.115E+00,&
-1.121E+00,1.129E+00,1.137E+00,1.147E+00,1.156E+00,1.164E+00,1.170E+00,1.175E+00,1.181E+00,1.188E+00,&
-1.196E+00,1.206E+00,1.216E+00,1.227E+00,1.239E+00,1.252E+00,1.267E+00,1.283E+00,1.297E+00,1.310E+00/
- data (k_Cw(i),i=901,1000)/&
-1.323E+00,1.336E+00,1.351E+00,1.370E+00,1.392E+00,1.416E+00,1.440E+00,1.465E+00,1.489E+00,1.511E+00,&
-1.532E+00,1.555E+00,1.580E+00,1.610E+00,1.642E+00,1.672E+00,1.701E+00,1.728E+00,1.758E+00,1.791E+00,&
-1.831E+00,1.872E+00,1.911E+00,1.943E+00,1.974E+00,2.007E+00,2.047E+00,2.098E+00,2.153E+00,2.203E+00,&
-2.243E+00,2.277E+00,2.312E+00,2.357E+00,2.415E+00,2.479E+00,2.540E+00,2.590E+00,2.631E+00,2.666E+00,&
-2.701E+00,2.738E+00,2.780E+00,2.829E+00,2.889E+00,2.960E+00,3.033E+00,3.097E+00,3.146E+00,3.181E+00,&
-3.226E+00,3.267E+00,3.319E+00,3.363E+00,3.412E+00,3.449E+00,3.504E+00,3.544E+00,3.600E+00,3.648E+00,&
-3.701E+00,3.752E+00,3.802E+00,3.871E+00,3.927E+00,3.985E+00,4.064E+00,4.125E+00,4.216E+00,4.302E+00,&
-4.389E+00,4.504E+00,4.630E+00,4.737E+00,4.904E+00,5.092E+00,5.260E+00,5.479E+00,5.720E+00,6.006E+00,&
-6.242E+00,6.580E+00,6.927E+00,7.313E+00,7.633E+00,8.089E+00,8.545E+00,9.030E+00,9.591E+00,1.002E+01,&
-1.063E+01,1.122E+01,1.184E+01,1.245E+01,1.316E+01,1.369E+01,1.434E+01,1.509E+01,1.578E+01,1.646E+01/
- data (k_Cw(i),i=1001,1100)/&
-1.714E+01,1.781E+01,1.854E+01,1.919E+01,1.980E+01,2.029E+01,2.089E+01,2.146E+01,2.202E+01,2.260E+01,&
-2.313E+01,2.360E+01,2.407E+01,2.450E+01,2.493E+01,2.533E+01,2.571E+01,2.606E+01,2.641E+01,2.673E+01,&
-2.701E+01,2.729E+01,2.756E+01,2.782E+01,2.806E+01,2.835E+01,2.856E+01,2.875E+01,2.892E+01,2.908E+01,&
-2.926E+01,2.940E+01,2.956E+01,2.966E+01,2.982E+01,2.993E+01,3.003E+01,3.014E+01,3.023E+01,3.029E+01,&
-3.036E+01,3.042E+01,3.046E+01,3.049E+01,3.052E+01,3.053E+01,3.055E+01,3.056E+01,3.056E+01,3.055E+01,&
-3.054E+01,3.051E+01,3.049E+01,3.045E+01,3.041E+01,3.035E+01,3.029E+01,3.023E+01,3.014E+01,3.006E+01,&
-2.998E+01,2.983E+01,2.971E+01,2.957E+01,2.936E+01,2.917E+01,2.899E+01,2.872E+01,2.851E+01,2.829E+01,&
-2.800E+01,2.777E+01,2.754E+01,2.722E+01,2.699E+01,2.664E+01,2.638E+01,2.611E+01,2.581E+01,2.555E+01,&
-2.522E+01,2.497E+01,2.468E+01,2.443E+01,2.413E+01,2.388E+01,2.364E+01,2.332E+01,2.307E+01,2.274E+01,&
-2.250E+01,2.218E+01,2.193E+01,2.163E+01,2.139E+01,2.107E+01,2.082E+01,2.052E+01,2.025E+01,2.001E+01/
- data (k_Cw(i),i=1101,1200)/&
-1.972E+01,1.951E+01,1.924E+01,1.900E+01,1.874E+01,1.847E+01,1.827E+01,1.802E+01,1.784E+01,1.758E+01,&
-1.734E+01,1.712E+01,1.688E+01,1.671E+01,1.647E+01,1.623E+01,1.606E+01,1.583E+01,1.562E+01,1.545E+01,&
-1.525E+01,1.504E+01,1.489E+01,1.468E+01,1.447E+01,1.432E+01,1.413E+01,1.395E+01,1.381E+01,1.364E+01,&
-1.348E+01,1.329E+01,1.316E+01,1.298E+01,1.282E+01,1.265E+01,1.254E+01,1.238E+01,1.223E+01,1.206E+01,&
-1.193E+01,1.181E+01,1.166E+01,1.152E+01,1.137E+01,1.126E+01,1.114E+01,1.100E+01,1.088E+01,1.075E+01,&
-1.064E+01,1.054E+01,1.044E+01,1.032E+01,1.022E+01,1.011E+01,1.001E+01,9.912E+00,9.839E+00,9.754E+00,&
-9.660E+00,9.563E+00,9.477E+00,9.383E+00,9.305E+00,9.202E+00,9.133E+00,9.047E+00,8.977E+00,8.898E+00,&
-8.820E+00,8.742E+00,8.665E+00,8.588E+00,8.509E+00,8.448E+00,8.364E+00,8.295E+00,8.234E+00,8.157E+00,&
-8.104E+00,8.036E+00,7.959E+00,7.890E+00,7.834E+00,7.773E+00,7.712E+00,7.654E+00,7.609E+00,7.548E+00,&
-7.495E+00,7.432E+00,7.374E+00,7.315E+00,7.252E+00,7.203E+00,7.164E+00,7.124E+00,7.084E+00,7.041E+00/
- data (k_Cw(i),i=1201,1300)/&
-6.987E+00,6.943E+00,6.910E+00,6.865E+00,6.828E+00,6.776E+00,6.742E+00,6.714E+00,6.695E+00,6.654E+00,&
-6.630E+00,6.599E+00,6.567E+00,6.526E+00,6.501E+00,6.474E+00,6.449E+00,6.420E+00,6.401E+00,6.363E+00,&
-6.345E+00,6.309E+00,6.282E+00,6.250E+00,6.214E+00,6.186E+00,6.163E+00,6.130E+00,6.121E+00,6.091E+00,&
-6.076E+00,6.053E+00,6.048E+00,6.016E+00,6.005E+00,5.982E+00,5.973E+00,5.947E+00,5.940E+00,5.919E+00,&
-5.911E+00,5.887E+00,5.875E+00,5.846E+00,5.826E+00,5.798E+00,5.787E+00,5.751E+00,5.746E+00,5.718E+00,&
-5.705E+00,5.685E+00,5.684E+00,5.657E+00,5.658E+00,5.644E+00,5.648E+00,5.626E+00,5.626E+00,5.619E+00,&
-5.618E+00,5.603E+00,5.614E+00,5.597E+00,5.603E+00,5.582E+00,5.584E+00,5.564E+00,5.563E+00,5.547E+00,&
-5.545E+00,5.536E+00,5.542E+00,5.529E+00,5.532E+00,5.525E+00,5.533E+00,5.528E+00,5.529E+00,5.516E+00,&
-5.524E+00,5.516E+00,5.526E+00,5.520E+00,5.520E+00,5.516E+00,5.522E+00,5.511E+00,5.527E+00,5.511E+00,&
-5.519E+00,5.515E+00,5.520E+00,5.510E+00,5.518E+00,5.523E+00,5.538E+00,5.535E+00,5.544E+00,5.557E+00/
- data (k_Cw(i),i=1301,1400)/&
-5.571E+00,5.583E+00,5.606E+00,5.607E+00,5.629E+00,5.636E+00,5.664E+00,5.670E+00,5.693E+00,5.702E+00,&
-5.733E+00,5.752E+00,5.766E+00,5.776E+00,5.797E+00,5.811E+00,5.829E+00,5.842E+00,5.877E+00,5.891E+00,&
-5.930E+00,5.945E+00,5.972E+00,5.999E+00,6.025E+00,6.051E+00,6.087E+00,6.096E+00,6.136E+00,6.166E+00,&
-6.198E+00,6.219E+00,6.256E+00,6.284E+00,6.335E+00,6.369E+00,6.392E+00,6.445E+00,6.493E+00,6.517E+00,&
-6.571E+00,6.617E+00,6.658E+00,6.689E+00,6.748E+00,6.796E+00,6.842E+00,6.897E+00,6.955E+00,7.003E+00,&
-7.054E+00,7.111E+00,7.179E+00,7.235E+00,7.274E+00,7.339E+00,7.414E+00,7.481E+00,7.536E+00,7.594E+00,&
-7.669E+00,7.734E+00,7.776E+00,7.833E+00,7.893E+00,7.952E+00,8.000E+00,8.045E+00,8.103E+00,8.155E+00,&
-8.205E+00,8.241E+00,8.264E+00,8.321E+00,8.352E+00,8.394E+00,8.430E+00,8.448E+00,8.477E+00,8.512E+00,&
-8.535E+00,8.562E+00,8.593E+00,8.618E+00,8.640E+00,8.670E+00,8.689E+00,8.720E+00,8.738E+00,8.755E+00,&
-8.777E+00,8.778E+00,8.778E+00,8.794E+00,8.805E+00,8.807E+00,8.809E+00,8.811E+00,8.799E+00,8.795E+00/
- data (k_Cw(i),i=1401,1500)/&
-8.789E+00,8.779E+00,8.767E+00,8.754E+00,8.750E+00,8.738E+00,8.739E+00,8.735E+00,8.744E+00,8.753E+00,&
-8.755E+00,8.780E+00,8.787E+00,8.790E+00,8.798E+00,8.794E+00,8.811E+00,8.820E+00,8.836E+00,8.845E+00,&
-8.854E+00,8.858E+00,8.868E+00,8.869E+00,8.884E+00,8.888E+00,8.900E+00,8.922E+00,8.951E+00,8.973E+00,&
-9.010E+00,9.034E+00,9.110E+00,9.146E+00,9.195E+00,9.259E+00,9.315E+00,9.380E+00,9.457E+00,9.535E+00,&
-9.633E+00,9.723E+00,9.824E+00,9.935E+00,1.005E+01,1.018E+01,1.031E+01,1.042E+01,1.059E+01,1.075E+01,&
-1.094E+01,1.110E+01,1.139E+01,1.160E+01,1.184E+01,1.208E+01,1.235E+01,1.266E+01,1.301E+01,1.334E+01,&
-1.375E+01,1.412E+01,1.478E+01,1.529E+01,1.586E+01,1.641E+01,1.709E+01,1.774E+01,1.853E+01,1.925E+01,&
-2.051E+01,2.148E+01,2.250E+01,2.367E+01,2.483E+01,2.622E+01,2.751E+01,2.972E+01,3.144E+01,3.317E+01,&
-3.504E+01,3.725E+01,3.927E+01,4.269E+01,4.530E+01,4.789E+01,5.060E+01,5.360E+01,5.761E+01,6.088E+01,&
-6.401E+01,6.720E+01,7.059E+01,7.497E+01,7.841E+01,8.157E+01,8.469E+01,8.752E+01,9.189E+01,9.456E+01/
- data (k_Cw(i),i=1501,1600)/&
-9.722E+01,9.995E+01,1.033E+02,1.057E+02,1.078E+02,1.097E+02,1.123E+02,1.140E+02,1.157E+02,1.171E+02,&
-1.190E+02,1.203E+02,1.215E+02,1.226E+02,1.240E+02,1.248E+02,1.257E+02,1.265E+02,1.275E+02,1.281E+02,&
-1.286E+02,1.292E+02,1.296E+02,1.299E+02,1.303E+02,1.304E+02,1.306E+02,1.306E+02,1.306E+02,1.305E+02,&
-1.304E+02,1.301E+02,1.299E+02,1.296E+02,1.291E+02,1.287E+02,1.282E+02,1.275E+02,1.270E+02,1.264E+02,&
-1.256E+02,1.249E+02,1.242E+02,1.232E+02,1.224E+02,1.216E+02,1.206E+02,1.197E+02,1.190E+02,1.178E+02,&
-1.170E+02,1.157E+02,1.149E+02,1.140E+02,1.128E+02,1.118E+02,1.110E+02,1.097E+02,1.088E+02,1.076E+02,&
-1.066E+02,1.058E+02,1.044E+02,1.036E+02,1.023E+02,1.014E+02,1.005E+02,9.928E+01,9.831E+01,9.711E+01,&
-9.631E+01,9.488E+01,9.412E+01,9.341E+01,9.206E+01,9.121E+01,9.009E+01,8.929E+01,8.804E+01,8.725E+01,&
-8.611E+01,8.532E+01,8.460E+01,8.336E+01,8.262E+01,8.151E+01,8.076E+01,7.973E+01,7.904E+01,7.800E+01,&
-7.723E+01,7.628E+01,7.557E+01,7.463E+01,7.392E+01,7.298E+01,7.234E+01,7.141E+01,7.082E+01,6.986E+01/
- data (k_Cw(i),i=1601,1700)/&
-6.924E+01,6.865E+01,6.779E+01,6.688E+01,6.634E+01,6.548E+01,6.490E+01,6.412E+01,6.358E+01,6.281E+01,&
-6.232E+01,6.156E+01,6.105E+01,6.029E+01,5.980E+01,5.907E+01,5.857E+01,5.788E+01,5.746E+01,5.680E+01,&
-5.632E+01,5.566E+01,5.503E+01,5.463E+01,5.398E+01,5.356E+01,5.294E+01,5.256E+01,5.197E+01,5.138E+01,&
-5.098E+01,5.047E+01,5.006E+01,4.951E+01,4.915E+01,4.860E+01,4.806E+01,4.773E+01,4.718E+01,4.685E+01,&
-4.635E+01,4.584E+01,4.553E+01,4.503E+01,4.473E+01,4.426E+01,4.379E+01,4.348E+01,4.302E+01,4.273E+01,&
-4.228E+01,4.186E+01,4.155E+01,4.115E+01,4.072E+01,4.046E+01,4.002E+01,3.975E+01,3.934E+01,3.894E+01,&
-3.868E+01,3.827E+01,3.790E+01,3.762E+01,3.726E+01,3.687E+01,3.664E+01,3.628E+01,3.592E+01,3.568E+01,&
-3.536E+01,3.499E+01,3.479E+01,3.447E+01,3.413E+01,3.392E+01,3.363E+01,3.329E+01,3.309E+01,3.280E+01,&
-3.251E+01,3.231E+01,3.200E+01,3.172E+01,3.154E+01,3.127E+01,3.097E+01,3.071E+01,3.052E+01,3.025E+01,&
-2.999E+01,2.980E+01,2.958E+01,2.929E+01,2.906E+01,2.888E+01,2.865E+01,2.840E+01,2.822E+01,2.799E+01/
- data (k_Cw(i),i=1701,1800)/&
-2.776E+01,2.751E+01,2.737E+01,2.713E+01,2.692E+01,2.671E+01,2.657E+01,2.637E+01,2.617E+01,2.597E+01,&
-2.584E+01,2.565E+01,2.547E+01,2.528E+01,2.518E+01,2.501E+01,2.485E+01,2.466E+01,2.455E+01,2.438E+01,&
-2.420E+01,2.401E+01,2.386E+01,2.374E+01,2.358E+01,2.342E+01,2.328E+01,2.318E+01,2.304E+01,2.289E+01,&
-2.275E+01,2.262E+01,2.252E+01,2.239E+01,2.227E+01,2.214E+01,2.202E+01,2.194E+01,2.183E+01,2.172E+01,&
-2.162E+01,2.150E+01,2.140E+01,2.133E+01,2.123E+01,2.112E+01,2.102E+01,2.091E+01,2.085E+01,2.077E+01,&
-2.067E+01,2.058E+01,2.049E+01,2.040E+01,2.033E+01,2.027E+01,2.018E+01,2.010E+01,2.002E+01,1.995E+01,&
-1.987E+01,1.983E+01,1.976E+01,1.969E+01,1.963E+01,1.958E+01,1.952E+01,1.946E+01,1.942E+01,1.938E+01,&
-1.933E+01,1.928E+01,1.923E+01,1.918E+01,1.913E+01,1.909E+01,1.904E+01,1.900E+01,1.896E+01,1.889E+01,&
-1.884E+01,1.880E+01,1.874E+01,1.870E+01,1.866E+01,1.860E+01,1.856E+01,1.852E+01,1.849E+01,1.845E+01,&
-1.843E+01,1.839E+01,1.837E+01,1.837E+01,1.838E+01,1.837E+01,1.836E+01,1.835E+01,1.835E+01,1.835E+01/
- data (k_Cw(i),i=1801,1900)/&
-1.834E+01,1.833E+01,1.833E+01,1.833E+01,1.834E+01,1.833E+01,1.834E+01,1.834E+01,1.832E+01,1.831E+01,&
-1.829E+01,1.828E+01,1.828E+01,1.829E+01,1.829E+01,1.830E+01,1.832E+01,1.835E+01,1.836E+01,1.837E+01,&
-1.837E+01,1.838E+01,1.837E+01,1.840E+01,1.842E+01,1.843E+01,1.846E+01,1.849E+01,1.851E+01,1.854E+01,&
-1.857E+01,1.858E+01,1.860E+01,1.864E+01,1.865E+01,1.867E+01,1.870E+01,1.875E+01,1.879E+01,1.883E+01,&
-1.888E+01,1.893E+01,1.897E+01,1.903E+01,1.908E+01,1.913E+01,1.919E+01,1.924E+01,1.929E+01,1.934E+01,&
-1.940E+01,1.946E+01,1.951E+01,1.961E+01,1.966E+01,1.973E+01,1.980E+01,1.985E+01,1.991E+01,1.997E+01,&
-2.004E+01,2.013E+01,2.020E+01,2.027E+01,2.038E+01,2.046E+01,2.054E+01,2.062E+01,2.071E+01,2.079E+01,&
-2.087E+01,2.094E+01,2.103E+01,2.111E+01,2.126E+01,2.135E+01,2.145E+01,2.153E+01,2.164E+01,2.174E+01,&
-2.183E+01,2.192E+01,2.202E+01,2.212E+01,2.223E+01,2.234E+01,2.245E+01,2.264E+01,2.276E+01,2.287E+01,&
-2.301E+01,2.313E+01,2.327E+01,2.338E+01,2.351E+01,2.363E+01,2.377E+01,2.385E+01,2.409E+01,2.418E+01/
- data (k_Cw(i),i=1901,2000)/&
-2.433E+01,2.444E+01,2.459E+01,2.470E+01,2.485E+01,2.495E+01,2.523E+01,2.534E+01,2.549E+01,2.564E+01,&
-2.579E+01,2.592E+01,2.611E+01,2.623E+01,2.653E+01,2.664E+01,2.681E+01,2.695E+01,2.712E+01,2.727E+01,&
-2.744E+01,2.756E+01,2.789E+01,2.802E+01,2.819E+01,2.838E+01,2.855E+01,2.869E+01,2.903E+01,2.916E+01,&
-2.934E+01,2.951E+01,2.969E+01,2.988E+01,3.002E+01,3.038E+01,3.054E+01,3.073E+01,3.092E+01,3.107E+01,&
-3.145E+01,3.160E+01,3.180E+01,3.199E+01,3.219E+01,3.238E+01,3.272E+01,3.290E+01,3.311E+01,3.331E+01,&
-3.348E+01,3.388E+01,3.404E+01,3.426E+01,3.447E+01,3.468E+01,3.486E+01,3.528E+01,3.547E+01,3.570E+01,&
-3.592E+01,3.610E+01,3.656E+01,3.675E+01,3.697E+01,3.716E+01,3.761E+01,3.777E+01,3.801E+01,3.827E+01,&
-3.848E+01,3.892E+01,3.913E+01,3.938E+01,3.957E+01,4.007E+01,4.024E+01,4.050E+01,4.078E+01,4.096E+01,&
-4.147E+01,4.167E+01,4.196E+01,4.220E+01,4.268E+01,4.293E+01,4.319E+01,4.342E+01,4.388E+01,4.413E+01,&
-4.439E+01,4.464E+01,4.513E+01,4.538E+01,4.565E+01,4.590E+01,4.638E+01,4.663E+01,4.690E+01,4.714E+01/
- data (k_Cw(i),i=2001,2101)/&
-4.764E+01,4.787E+01,4.817E+01,4.836E+01,4.894E+01,4.917E+01,4.940E+01,4.997E+01,5.021E+01,5.050E+01,&
-5.078E+01,5.131E+01,5.159E+01,5.191E+01,5.246E+01,5.275E+01,5.310E+01,5.339E+01,5.401E+01,5.429E+01,&
-5.462E+01,5.523E+01,5.555E+01,5.588E+01,5.618E+01,5.686E+01,5.713E+01,5.746E+01,5.808E+01,5.841E+01,&
-5.878E+01,5.946E+01,5.976E+01,6.009E+01,6.039E+01,6.107E+01,6.146E+01,6.180E+01,6.250E+01,6.286E+01,&
-6.316E+01,6.392E+01,6.427E+01,6.462E+01,6.539E+01,6.574E+01,6.609E+01,6.685E+01,6.725E+01,6.759E+01,&
-6.842E+01,6.881E+01,6.918E+01,7.009E+01,7.046E+01,7.084E+01,7.174E+01,7.212E+01,7.249E+01,7.340E+01,&
-7.381E+01,7.423E+01,7.514E+01,7.554E+01,7.591E+01,7.675E+01,7.718E+01,7.801E+01,7.839E+01,7.880E+01,&
-7.970E+01,8.005E+01,8.046E+01,8.136E+01,8.173E+01,8.215E+01,8.295E+01,8.338E+01,8.423E+01,8.474E+01,&
-8.510E+01,8.584E+01,8.617E+01,8.637E+01,8.708E+01,8.762E+01,8.884E+01,8.904E+01,8.945E+01,9.055E+01,&
-9.089E+01,9.134E+01,9.204E+01,9.244E+01,9.308E+01,9.333E+01,9.358E+01,9.416E+01,9.448E+01,9.499E+01,&
-9.530E+01/
-
-
-! ********************************************************************************
-! Specific absorption coefficient of dry matter
-! ********************************************************************************
-
- data (k_Cm(i),i=1,100)/&
-1.097E+02,1.037E+02,9.798E+01,9.244E+01,8.713E+01,8.231E+01,7.806E+01,7.404E+01,7.013E+01,6.654E+01,&
-6.300E+01,5.954E+01,5.616E+01,5.301E+01,5.001E+01,4.723E+01,4.463E+01,4.220E+01,3.996E+01,3.780E+01,&
-3.567E+01,3.362E+01,3.170E+01,2.993E+01,2.832E+01,2.679E+01,2.535E+01,2.402E+01,2.276E+01,2.150E+01,&
-2.024E+01,1.901E+01,1.785E+01,1.676E+01,1.575E+01,1.481E+01,1.392E+01,1.312E+01,1.233E+01,1.161E+01,&
-1.096E+01,1.041E+01,9.924E+00,9.410E+00,8.947E+00,8.508E+00,8.087E+00,7.640E+00,7.269E+00,6.939E+00,&
-6.660E+00,6.422E+00,6.222E+00,6.010E+00,5.782E+00,5.573E+00,5.370E+00,5.173E+00,4.946E+00,4.761E+00,&
-4.575E+00,4.419E+00,4.259E+00,4.117E+00,4.006E+00,3.945E+00,3.853E+00,3.784E+00,3.671E+00,3.554E+00,&
-3.462E+00,3.364E+00,3.282E+00,3.184E+00,3.102E+00,3.051E+00,2.983E+00,2.947E+00,2.913E+00,2.869E+00,&
-2.803E+00,2.777E+00,2.751E+00,2.726E+00,2.702E+00,2.679E+00,2.656E+00,2.634E+00,2.613E+00,2.593E+00,&
-2.573E+00,2.554E+00,2.536E+00,2.519E+00,2.502E+00,2.486E+00,2.471E+00,2.457E+00,2.443E+00,2.430E+00/
- data (k_Cm(i),i=101,200)/&
-2.417E+00,2.405E+00,2.394E+00,2.384E+00,2.374E+00,2.365E+00,2.356E+00,2.348E+00,2.341E+00,2.334E+00,&
-2.328E+00,2.323E+00,2.318E+00,2.314E+00,2.310E+00,2.307E+00,2.304E+00,2.303E+00,2.301E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00/
- data (k_Cm(i),i=201,300)/&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00/
- data (k_Cm(i),i=301,400)/&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00/
- data (k_Cm(i),i=401,500)/&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00/
- data (k_Cm(i),i=501,600)/&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00/
- data (k_Cm(i),i=601,700)/&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00/
- data (k_Cm(i),i=701,800)/&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,&
-2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00,2.300E+00/
- data (k_Cm(i),i=801,900)/&
-2.300E+00,2.300E+00,2.301E+00,2.301E+00,2.302E+00,2.303E+00,2.305E+00,2.306E+00,2.308E+00,2.310E+00,&
-2.312E+00,2.315E+00,2.318E+00,2.320E+00,2.323E+00,2.327E+00,2.330E+00,2.334E+00,2.338E+00,2.342E+00,&
-2.346E+00,2.350E+00,2.354E+00,2.359E+00,2.364E+00,2.369E+00,2.374E+00,2.379E+00,2.384E+00,2.389E+00,&
-2.395E+00,2.400E+00,2.406E+00,2.412E+00,2.418E+00,2.424E+00,2.430E+00,2.436E+00,2.442E+00,2.448E+00,&
-2.454E+00,2.461E+00,2.467E+00,2.473E+00,2.480E+00,2.486E+00,2.493E+00,2.499E+00,2.506E+00,2.512E+00,&
-2.519E+00,2.525E+00,2.532E+00,2.538E+00,2.545E+00,2.551E+00,2.558E+00,2.564E+00,2.570E+00,2.577E+00,&
-2.583E+00,2.589E+00,2.595E+00,2.601E+00,2.607E+00,2.613E+00,2.619E+00,2.625E+00,2.631E+00,2.636E+00,&
-2.642E+00,2.647E+00,2.652E+00,2.657E+00,2.662E+00,2.667E+00,2.672E+00,2.677E+00,2.681E+00,2.685E+00,&
-2.689E+00,2.693E+00,2.697E+00,2.701E+00,2.704E+00,2.708E+00,2.711E+00,2.714E+00,2.716E+00,2.719E+00,&
-2.721E+00,2.723E+00,2.725E+00,2.726E+00,2.728E+00,2.729E+00,2.730E+00,2.730E+00,2.731E+00,2.731E+00/
- data (k_Cm(i),i=901,1000)/&
-2.721E+00,2.736E+00,2.728E+00,2.719E+00,2.712E+00,2.698E+00,2.702E+00,2.691E+00,2.688E+00,2.686E+00,&
-2.682E+00,2.676E+00,2.675E+00,2.677E+00,2.670E+00,2.671E+00,2.668E+00,2.670E+00,2.674E+00,2.674E+00,&
-2.683E+00,2.674E+00,2.675E+00,2.682E+00,2.683E+00,2.665E+00,2.661E+00,2.668E+00,2.670E+00,2.674E+00,&
-2.664E+00,2.661E+00,2.665E+00,2.668E+00,2.681E+00,2.686E+00,2.684E+00,2.697E+00,2.712E+00,2.720E+00,&
-2.717E+00,2.726E+00,2.744E+00,2.743E+00,2.751E+00,2.763E+00,2.778E+00,2.793E+00,2.818E+00,2.835E+00,&
-2.865E+00,2.879E+00,2.899E+00,2.918E+00,2.936E+00,2.953E+00,2.966E+00,2.977E+00,2.981E+00,2.942E+00,&
-2.888E+00,2.864E+00,2.877E+00,2.886E+00,2.888E+00,2.891E+00,2.899E+00,2.887E+00,2.884E+00,2.900E+00,&
-2.929E+00,2.969E+00,3.014E+00,3.053E+00,3.075E+00,3.111E+00,3.128E+00,3.130E+00,3.103E+00,3.051E+00,&
-2.980E+00,2.941E+00,2.920E+00,2.931E+00,2.950E+00,2.979E+00,3.025E+00,3.048E+00,3.066E+00,3.087E+00,&
-3.099E+00,3.090E+00,3.088E+00,3.081E+00,3.086E+00,3.071E+00,3.065E+00,3.069E+00,3.067E+00,3.085E+00/
- data (k_Cm(i),i=1001,1100)/&
-3.094E+00,3.110E+00,3.136E+00,3.149E+00,3.158E+00,3.191E+00,3.230E+00,3.266E+00,3.298E+00,3.356E+00,&
-3.419E+00,3.476E+00,3.534E+00,3.584E+00,3.632E+00,3.708E+00,3.775E+00,3.847E+00,3.931E+00,3.987E+00,&
-4.071E+00,4.156E+00,4.242E+00,4.320E+00,4.395E+00,4.480E+00,4.561E+00,4.638E+00,4.708E+00,4.782E+00,&
-4.846E+00,4.906E+00,4.974E+00,5.026E+00,5.071E+00,5.131E+00,5.179E+00,5.220E+00,5.271E+00,5.322E+00,&
-5.358E+00,5.403E+00,5.441E+00,5.460E+00,5.481E+00,5.500E+00,5.523E+00,5.548E+00,5.560E+00,5.575E+00,&
-5.582E+00,5.597E+00,5.611E+00,5.639E+00,5.653E+00,5.675E+00,5.682E+00,5.685E+00,5.680E+00,5.689E+00,&
-5.711E+00,5.723E+00,5.715E+00,5.716E+00,5.732E+00,5.741E+00,5.743E+00,5.752E+00,5.745E+00,5.744E+00,&
-5.757E+00,5.766E+00,5.781E+00,5.787E+00,5.798E+00,5.810E+00,5.808E+00,5.815E+00,5.825E+00,5.824E+00,&
-5.827E+00,5.854E+00,5.878E+00,5.900E+00,5.908E+00,5.922E+00,5.940E+00,5.962E+00,5.963E+00,5.966E+00,&
-5.982E+00,5.990E+00,5.994E+00,6.016E+00,6.014E+00,6.025E+00,6.008E+00,6.022E+00,6.021E+00,6.027E+00/
- data (k_Cm(i),i=1101,1200)/&
-6.027E+00,6.035E+00,6.025E+00,6.009E+00,5.990E+00,5.987E+00,5.984E+00,5.971E+00,5.971E+00,5.973E+00,&
-5.951E+00,5.952E+00,5.939E+00,5.933E+00,5.931E+00,5.925E+00,5.909E+00,5.897E+00,5.884E+00,5.881E+00,&
-5.876E+00,5.856E+00,5.843E+00,5.830E+00,5.818E+00,5.807E+00,5.799E+00,5.792E+00,5.776E+00,5.783E+00,&
-5.776E+00,5.767E+00,5.762E+00,5.769E+00,5.756E+00,5.762E+00,5.737E+00,5.740E+00,5.757E+00,5.756E+00,&
-5.751E+00,5.754E+00,5.751E+00,5.750E+00,5.744E+00,5.744E+00,5.754E+00,5.744E+00,5.735E+00,5.740E+00,&
-5.732E+00,5.728E+00,5.731E+00,5.724E+00,5.724E+00,5.715E+00,5.697E+00,5.693E+00,5.703E+00,5.700E+00,&
-5.713E+00,5.728E+00,5.731E+00,5.735E+00,5.743E+00,5.754E+00,5.753E+00,5.744E+00,5.746E+00,5.746E+00,&
-5.746E+00,5.750E+00,5.749E+00,5.743E+00,5.740E+00,5.747E+00,5.739E+00,5.755E+00,5.749E+00,5.753E+00,&
-5.745E+00,5.732E+00,5.735E+00,5.724E+00,5.725E+00,5.712E+00,5.702E+00,5.700E+00,5.700E+00,5.686E+00,&
-5.685E+00,5.672E+00,5.659E+00,5.627E+00,5.613E+00,5.590E+00,5.591E+00,5.563E+00,5.552E+00,5.525E+00/
- data (k_Cm(i),i=1201,1300)/&
-5.517E+00,5.506E+00,5.494E+00,5.459E+00,5.450E+00,5.438E+00,5.428E+00,5.407E+00,5.391E+00,5.382E+00,&
-5.376E+00,5.358E+00,5.347E+00,5.333E+00,5.319E+00,5.301E+00,5.292E+00,5.284E+00,5.274E+00,5.258E+00,&
-5.253E+00,5.249E+00,5.233E+00,5.222E+00,5.211E+00,5.204E+00,5.210E+00,5.200E+00,5.193E+00,5.186E+00,&
-5.177E+00,5.175E+00,5.178E+00,5.173E+00,5.152E+00,5.135E+00,5.139E+00,5.128E+00,5.121E+00,5.114E+00,&
-5.126E+00,5.107E+00,5.104E+00,5.100E+00,5.109E+00,5.114E+00,5.112E+00,5.128E+00,5.137E+00,5.131E+00,&
-5.152E+00,5.175E+00,5.194E+00,5.200E+00,5.244E+00,5.257E+00,5.273E+00,5.289E+00,5.335E+00,5.366E+00,&
-5.389E+00,5.427E+00,5.453E+00,5.490E+00,5.522E+00,5.562E+00,5.605E+00,5.652E+00,5.698E+00,5.744E+00,&
-5.743E+00,5.642E+00,5.689E+00,5.722E+00,5.755E+00,5.798E+00,5.848E+00,5.875E+00,5.918E+00,5.971E+00,&
-6.022E+00,6.061E+00,6.116E+00,6.173E+00,6.214E+00,6.266E+00,6.319E+00,6.382E+00,6.426E+00,6.486E+00,&
-6.542E+00,6.579E+00,6.616E+00,6.666E+00,6.728E+00,6.771E+00,6.807E+00,6.858E+00,6.908E+00,6.959E+00/
- data (k_Cm(i),i=1301,1400)/&
-7.006E+00,7.052E+00,7.093E+00,7.136E+00,7.164E+00,7.199E+00,7.232E+00,7.266E+00,7.315E+00,7.340E+00,&
-7.361E+00,7.399E+00,7.440E+00,7.473E+00,7.505E+00,7.525E+00,7.557E+00,7.579E+00,7.604E+00,7.633E+00,&
-7.653E+00,7.674E+00,7.691E+00,7.699E+00,7.708E+00,7.721E+00,7.730E+00,7.729E+00,7.712E+00,7.702E+00,&
-7.695E+00,7.670E+00,7.644E+00,7.618E+00,7.585E+00,7.555E+00,7.517E+00,7.479E+00,7.451E+00,7.435E+00,&
-7.408E+00,7.379E+00,7.363E+00,7.347E+00,7.332E+00,7.332E+00,7.332E+00,7.304E+00,7.295E+00,7.296E+00,&
-7.291E+00,7.292E+00,7.292E+00,7.281E+00,7.283E+00,7.264E+00,7.250E+00,7.240E+00,7.228E+00,7.210E+00,&
-7.186E+00,7.164E+00,7.143E+00,7.114E+00,7.101E+00,7.069E+00,7.038E+00,7.003E+00,6.974E+00,6.928E+00,&
-6.889E+00,6.839E+00,6.793E+00,6.764E+00,6.729E+00,6.694E+00,6.662E+00,6.613E+00,6.572E+00,6.546E+00,&
-6.522E+00,6.507E+00,6.482E+00,6.484E+00,6.479E+00,6.494E+00,6.496E+00,6.491E+00,6.461E+00,6.440E+00,&
-6.430E+00,6.413E+00,6.421E+00,6.399E+00,6.379E+00,6.365E+00,6.372E+00,6.346E+00,6.321E+00,6.310E+00/
- data (k_Cm(i),i=1401,1500)/&
-6.314E+00,6.282E+00,6.277E+00,6.270E+00,6.258E+00,6.242E+00,6.234E+00,6.221E+00,6.231E+00,6.221E+00,&
-6.205E+00,6.193E+00,6.192E+00,6.179E+00,6.159E+00,6.143E+00,6.143E+00,6.120E+00,6.098E+00,6.087E+00,&
-6.063E+00,6.056E+00,6.053E+00,6.040E+00,6.044E+00,6.007E+00,5.996E+00,5.994E+00,5.997E+00,5.975E+00,&
-5.954E+00,5.946E+00,5.927E+00,5.914E+00,5.890E+00,5.873E+00,5.832E+00,5.794E+00,5.768E+00,5.728E+00,&
-5.681E+00,5.680E+00,5.655E+00,5.648E+00,5.620E+00,5.594E+00,5.567E+00,5.557E+00,5.552E+00,5.553E+00,&
-5.539E+00,5.524E+00,5.507E+00,5.505E+00,5.487E+00,5.474E+00,5.462E+00,5.450E+00,5.448E+00,5.441E+00,&
-5.440E+00,5.442E+00,5.450E+00,5.466E+00,5.461E+00,5.452E+00,5.445E+00,5.412E+00,5.379E+00,5.280E+00,&
-5.228E+00,5.199E+00,5.171E+00,5.139E+00,5.124E+00,5.112E+00,5.129E+00,5.145E+00,5.173E+00,5.176E+00,&
-5.185E+00,5.183E+00,5.200E+00,5.202E+00,5.204E+00,5.224E+00,5.244E+00,5.297E+00,5.318E+00,5.376E+00,&
-5.441E+00,5.491E+00,5.561E+00,5.629E+00,5.687E+00,5.757E+00,5.830E+00,5.911E+00,5.998E+00,6.090E+00/
- data (k_Cm(i),i=1501,1600)/&
-6.193E+00,6.330E+00,6.449E+00,6.581E+00,6.710E+00,6.838E+00,6.970E+00,7.107E+00,7.238E+00,7.388E+00,&
-7.506E+00,7.635E+00,7.767E+00,7.874E+00,7.977E+00,8.071E+00,8.150E+00,8.220E+00,8.292E+00,8.350E+00,&
-8.449E+00,8.521E+00,8.583E+00,8.666E+00,8.723E+00,8.759E+00,8.821E+00,8.864E+00,8.909E+00,8.941E+00,&
-8.949E+00,8.955E+00,8.983E+00,9.022E+00,9.043E+00,9.044E+00,9.028E+00,9.034E+00,9.052E+00,9.048E+00,&
-9.041E+00,9.037E+00,9.036E+00,9.035E+00,9.021E+00,9.016E+00,9.008E+00,8.970E+00,8.974E+00,8.953E+00,&
-8.957E+00,8.937E+00,8.923E+00,8.912E+00,8.895E+00,8.891E+00,8.880E+00,8.867E+00,8.855E+00,8.852E+00,&
-8.861E+00,8.864E+00,8.876E+00,8.869E+00,8.873E+00,8.855E+00,8.828E+00,8.839E+00,8.855E+00,8.856E+00,&
-8.833E+00,8.842E+00,8.844E+00,8.830E+00,8.808E+00,8.818E+00,8.807E+00,8.797E+00,8.794E+00,8.791E+00,&
-8.795E+00,8.772E+00,8.754E+00,8.759E+00,8.760E+00,8.746E+00,8.762E+00,8.778E+00,8.790E+00,8.795E+00,&
-8.811E+00,8.848E+00,8.874E+00,8.885E+00,8.913E+00,8.944E+00,8.981E+00,8.988E+00,9.001E+00,9.034E+00/
- data (k_Cm(i),i=1601,1700)/&
-9.076E+00,9.111E+00,9.141E+00,9.171E+00,9.214E+00,9.255E+00,9.304E+00,9.356E+00,9.406E+00,9.448E+00,&
-9.516E+00,9.578E+00,9.638E+00,9.692E+00,9.763E+00,9.845E+00,9.953E+00,1.004E+01,1.015E+01,1.027E+01,&
-1.039E+01,1.052E+01,1.063E+01,1.077E+01,1.091E+01,1.103E+01,1.119E+01,1.135E+01,1.150E+01,1.166E+01,&
-1.181E+01,1.201E+01,1.217E+01,1.235E+01,1.251E+01,1.269E+01,1.287E+01,1.307E+01,1.325E+01,1.346E+01,&
-1.364E+01,1.384E+01,1.404E+01,1.423E+01,1.441E+01,1.461E+01,1.481E+01,1.499E+01,1.518E+01,1.534E+01,&
-1.554E+01,1.571E+01,1.591E+01,1.607E+01,1.622E+01,1.637E+01,1.653E+01,1.667E+01,1.678E+01,1.690E+01,&
-1.698E+01,1.709E+01,1.718E+01,1.725E+01,1.734E+01,1.739E+01,1.748E+01,1.755E+01,1.761E+01,1.767E+01,&
-1.771E+01,1.777E+01,1.783E+01,1.787E+01,1.794E+01,1.795E+01,1.799E+01,1.805E+01,1.809E+01,1.813E+01,&
-1.820E+01,1.827E+01,1.830E+01,1.835E+01,1.841E+01,1.846E+01,1.852E+01,1.856E+01,1.861E+01,1.866E+01,&
-1.871E+01,1.876E+01,1.881E+01,1.885E+01,1.890E+01,1.896E+01,1.903E+01,1.907E+01,1.911E+01,1.916E+01/
- data (k_Cm(i),i=1701,1800)/&
-1.921E+01,1.927E+01,1.929E+01,1.932E+01,1.935E+01,1.936E+01,1.940E+01,1.943E+01,1.947E+01,1.951E+01,&
-1.953E+01,1.956E+01,1.960E+01,1.961E+01,1.962E+01,1.966E+01,1.966E+01,1.969E+01,1.970E+01,1.972E+01,&
-1.974E+01,1.974E+01,1.976E+01,1.978E+01,1.979E+01,1.982E+01,1.982E+01,1.983E+01,1.986E+01,1.988E+01,&
-1.989E+01,1.989E+01,1.993E+01,1.993E+01,1.997E+01,2.001E+01,2.001E+01,2.004E+01,2.007E+01,2.011E+01,&
-2.014E+01,2.015E+01,2.017E+01,2.019E+01,2.022E+01,2.023E+01,2.023E+01,2.026E+01,2.027E+01,2.029E+01,&
-2.028E+01,2.027E+01,2.029E+01,2.029E+01,2.031E+01,2.029E+01,2.026E+01,2.030E+01,2.031E+01,2.030E+01,&
-2.030E+01,2.031E+01,2.029E+01,2.029E+01,2.026E+01,2.025E+01,2.023E+01,2.020E+01,2.016E+01,2.015E+01,&
-2.012E+01,2.009E+01,2.007E+01,2.003E+01,1.998E+01,1.996E+01,1.991E+01,1.986E+01,1.980E+01,1.975E+01,&
-1.969E+01,1.964E+01,1.959E+01,1.953E+01,1.947E+01,1.941E+01,1.936E+01,1.931E+01,1.922E+01,1.918E+01,&
-1.912E+01,1.906E+01,1.899E+01,1.890E+01,1.885E+01,1.881E+01,1.875E+01,1.867E+01,1.861E+01,1.858E+01/
- data (k_Cm(i),i=1801,1900)/&
-1.852E+01,1.846E+01,1.840E+01,1.834E+01,1.829E+01,1.824E+01,1.819E+01,1.813E+01,1.807E+01,1.803E+01,&
-1.798E+01,1.792E+01,1.788E+01,1.782E+01,1.780E+01,1.775E+01,1.773E+01,1.768E+01,1.768E+01,1.766E+01,&
-1.763E+01,1.762E+01,1.763E+01,1.764E+01,1.764E+01,1.766E+01,1.770E+01,1.774E+01,1.779E+01,1.786E+01,&
-1.795E+01,1.803E+01,1.814E+01,1.825E+01,1.836E+01,1.851E+01,1.866E+01,1.881E+01,1.895E+01,1.913E+01,&
-1.932E+01,1.951E+01,1.972E+01,1.994E+01,2.017E+01,2.040E+01,2.065E+01,2.089E+01,2.114E+01,2.140E+01,&
-2.163E+01,2.186E+01,2.210E+01,2.237E+01,2.262E+01,2.290E+01,2.313E+01,2.339E+01,2.361E+01,2.387E+01,&
-2.410E+01,2.435E+01,2.455E+01,2.479E+01,2.499E+01,2.521E+01,2.541E+01,2.562E+01,2.583E+01,2.605E+01,&
-2.626E+01,2.643E+01,2.657E+01,2.674E+01,2.689E+01,2.701E+01,2.718E+01,2.732E+01,2.742E+01,2.754E+01,&
-2.763E+01,2.777E+01,2.794E+01,2.804E+01,2.821E+01,2.836E+01,2.850E+01,2.863E+01,2.878E+01,2.896E+01,&
-2.913E+01,2.922E+01,2.937E+01,2.947E+01,2.960E+01,2.970E+01,2.982E+01,2.997E+01,3.007E+01,3.018E+01/
- data (k_Cm(i),i=1901,2000)/&
-3.028E+01,3.040E+01,3.053E+01,3.060E+01,3.066E+01,3.070E+01,3.076E+01,3.078E+01,3.075E+01,3.074E+01,&
-3.072E+01,3.065E+01,3.058E+01,3.051E+01,3.045E+01,3.034E+01,3.029E+01,3.023E+01,3.015E+01,3.004E+01,&
-3.000E+01,2.998E+01,2.991E+01,2.986E+01,2.984E+01,2.981E+01,2.976E+01,2.973E+01,2.976E+01,2.977E+01,&
-2.979E+01,2.981E+01,2.985E+01,2.989E+01,2.999E+01,3.000E+01,3.005E+01,3.007E+01,3.011E+01,3.017E+01,&
-3.022E+01,3.023E+01,3.029E+01,3.028E+01,3.029E+01,3.027E+01,3.027E+01,3.024E+01,3.019E+01,3.010E+01,&
-3.010E+01,3.003E+01,2.993E+01,2.983E+01,2.984E+01,2.975E+01,2.966E+01,2.958E+01,2.948E+01,2.930E+01,&
-2.926E+01,2.920E+01,2.913E+01,2.902E+01,2.890E+01,2.882E+01,2.873E+01,2.873E+01,2.870E+01,2.865E+01,&
-2.858E+01,2.854E+01,2.851E+01,2.846E+01,2.838E+01,2.834E+01,2.823E+01,2.820E+01,2.817E+01,2.809E+01,&
-2.803E+01,2.803E+01,2.801E+01,2.794E+01,2.791E+01,2.790E+01,2.784E+01,2.779E+01,2.781E+01,2.782E+01,&
-2.781E+01,2.781E+01,2.783E+01,2.785E+01,2.785E+01,2.782E+01,2.782E+01,2.780E+01,2.780E+01,2.778E+01/
- data (k_Cm(i),i=2001,2101)/&
-2.783E+01,2.784E+01,2.788E+01,2.791E+01,2.799E+01,2.804E+01,2.812E+01,2.811E+01,2.819E+01,2.819E+01,&
-2.817E+01,2.817E+01,2.831E+01,2.837E+01,2.848E+01,2.853E+01,2.859E+01,2.870E+01,2.874E+01,2.887E+01,&
-2.898E+01,2.910E+01,2.923E+01,2.934E+01,2.944E+01,2.959E+01,2.973E+01,2.987E+01,3.002E+01,3.016E+01,&
-3.035E+01,3.043E+01,3.064E+01,3.084E+01,3.098E+01,3.122E+01,3.132E+01,3.152E+01,3.165E+01,3.184E+01,&
-3.204E+01,3.221E+01,3.233E+01,3.255E+01,3.282E+01,3.315E+01,3.339E+01,3.360E+01,3.384E+01,3.410E+01,&
-3.426E+01,3.452E+01,3.473E+01,3.497E+01,3.519E+01,3.539E+01,3.561E+01,3.579E+01,3.604E+01,3.618E+01,&
-3.636E+01,3.660E+01,3.675E+01,3.682E+01,3.699E+01,3.711E+01,3.735E+01,3.758E+01,3.792E+01,3.796E+01,&
-3.812E+01,3.822E+01,3.833E+01,3.856E+01,3.874E+01,3.877E+01,3.884E+01,3.883E+01,3.879E+01,3.886E+01,&
-3.900E+01,3.900E+01,3.906E+01,3.905E+01,3.916E+01,3.928E+01,3.948E+01,3.943E+01,3.951E+01,3.964E+01,&
-3.953E+01,3.960E+01,3.958E+01,3.954E+01,3.940E+01,3.936E+01,3.917E+01,3.926E+01,3.893E+01,3.921E+01,&
-3.871E+01/
-
-end
\ No newline at end of file
diff --git a/src/models/pktestProspect.cc b/src/models/pktestProspect.cc
deleted file mode 100644
index 00ff3d6..0000000
--- a/src/models/pktestProspect.cc
+++ /dev/null
@@ -1,76 +0,0 @@
-/**********************************************************************
-pktestProspect: example program how to use class Prospect
-Copyright (C) 2008-2013 Pieter Kempeneers
-
-This file is part of pktools
-
-pktools is free software: you can redistribute it and/or modify
-it under the terms of the GNU General Public License as published by
-the Free Software Foundation, either version 3 of the License, or
-(at your option) any later version.
-
-pktools is distributed in the hope that it will be useful,
-but WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-GNU General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with pktools. If not, see <http://www.gnu.org/licenses/>.
-***********************************************************************/
-#include <iostream>
-#include <string>
-#include "base/Optionpk.h"
-#include "algorithms/Filter.h"
-#include "fileclasses/FileReaderAscii.h"
-#include "Prospect.h"
-
-int main(int argc, char *argv[])
-{
- Optionpk<double> N_opt("N","N","parameter N for Prospect model",1.2);
- Optionpk<double> Cab_opt("Cab","Cab","parameter Cab for Prospect model (g.cm-2)",30);
- Optionpk<double> Car_opt("Car","Car","parameter Car for Prospect model (g.cm-2)",10);
- Optionpk<double> Cbrown_opt("Cbrown","Cbrown","parameter Cbrown for Prospect model (arbitrary units)",0);
- Optionpk<double> Cw_opt("Cw","Cw","parameter Cw for Prospect model (cm)",0.015);
- Optionpk<double> Cm_opt("Cm","Cm","parameter Cm for Prospect model (g.cm-2)",0.009);
- Optionpk<double> step_opt("s","step","step size for interpolation",0.01);
- Optionpk<double> fwhm_opt("fwhm", "fwhm", "list of full width half to apply spectral filtering (-fwhm band1 -fwhm band2 ...)");
- Optionpk<double> wavelengthIn_opt("win", "wavelengthIn", "list of wavelengths in input spectrum (-w band1 -w band2 ...)");
- Optionpk<double> wavelengthOut_opt("wout", "wavelengthOut", "list of wavelengths in output spectrum (-w band1 -w band2 ...)");
- Optionpk<short> verbose_opt("v","verbose","verbose mode",0);
-
- bool doProcess;//stop process when program was invoked with help option (-h --help)
- try{
- doProcess=N_opt.retrieveOption(argc,argv);
- Cab_opt.retrieveOption(argc,argv);
- Car_opt.retrieveOption(argc,argv);
- Cbrown_opt.retrieveOption(argc,argv);
- Cw_opt.retrieveOption(argc,argv);
- Cm_opt.retrieveOption(argc,argv);
- fwhm_opt.retrieveOption(argc,argv);
- wavelengthIn_opt.retrieveOption(argc,argv);
- wavelengthOut_opt.retrieveOption(argc,argv);
- verbose_opt.retrieveOption(argc,argv);
- }
- catch(std::string predefinedString){//command line option contained license or version
- std::cout << predefinedString << std::endl;//report the predefined string to stdout
- exit(0);//stop processing
- }
- if(!doProcess){//command line option contained help option
- std::cout << "short option -h shows basic options only, use long option --help to show all options" << std::endl;//provide extra details for help to the user
- exit(0);//stop processing
- }
-
- prospect::Prospect prospectModel;
- prospectModel.setN(N_opt[0]);
- prospectModel.setCab(Cab_opt[0]);
- prospectModel.setCar(Car_opt[0]);
- prospectModel.setCbrown(Cbrown_opt[0]);
- prospectModel.setCw(Cw_opt[0]);
- prospectModel.setCm(Cm_opt[0]);
- vector<double> leafReflectance;
- vector<double> leafTransmittance;
- prospectModel.getLeafSpectrum(leafReflectance,leafTransmittance);
- assert(leafReflectance.size()==leafTransmittance.size());
- for(int index=0;index<leafReflectance.size();++index)
- std::cout << 400+index << " " << leafReflectance[index] << " " << leafTransmittance[index] << std::endl;
-}
diff --git a/src/models/prospect_5B.f90 b/src/models/prospect_5B.f90
deleted file mode 100644
index 331a593..0000000
--- a/src/models/prospect_5B.f90
+++ /dev/null
@@ -1,132 +0,0 @@
-! ********************************************************************************
-! prospect_5B.f90
-! ********************************************************************************
-!
-! Jean-Baptiste FERET
-!
-! Department of Global Ecology / Carnegie Institution for Sciences
-! 260 Panama Street
-! Stanford, CA 94305, USA
-! E-mail: jbferet at stanford.edu
-!
-! St�phane JACQUEMOUD
-!
-! Universit� Paris Diderot / Institut de Physique du Globe de Paris
-! 35 rue H�l�ne Brion
-! 75013 Paris, France
-! E-mail: jacquemoud at ipgp.fr
-!
-! http://teledetection.ipgp.fr/prosail/
-!
-! ********************************************************************************
-! F�ret J.B., Fran�ois C., Asner G.P., Gitelson A.A., Martin R.E., Bidel L.P.R.,
-! Ustin S.L., le Maire G., Jacquemoud S. (2008), PROSPECT-4 and 5: Advances in the
-! leaf optical properties model separating photosynthetic pigments, Remote Sensing
-! of Environment, 112:3030-3043.
-! Jacquemoud S., Ustin S.L., Verdebout J., Schmuck G., Andreoli G., Hosgood B.
-! (1996), Estimating leaf biochemistry using the PROSPECT leaf optical properties
-! model, Remote Sensing of Environment, 56:194-202.
-! Jacquemoud S., Baret F. (1990), PROSPECT: a model of leaf optical properties
-! spectra, Remote Sensing of Environment, 34:75-91.
-! ********************************************************************************
-! version 5.02 (25 July 2011)
-! ********************************************************************************
-
-subroutine prospect_5B(N,Cab,Car,Cbrown,Cw,Cm,RT)
-
-use dataSpec_P5B
-implicit none
-
-real(8), intent(in) :: N,Cab,Car,Cbrown,Cw,Cm
-real(8), intent(out) :: RT(nw,2)
-
-real(8) :: k(nw), tau(nw), xx(nw), yy(nw)
-real(8) :: x1(nw),x2(nw),x3(nw),x4(nw),x5(nw),x6(nw)
-real(8) :: theta1, theta2, t1(nw),t2(nw)
-real(8) :: r(nw),t(nw),ra(nw),ta(nw)
-real(8) :: delta(nw),beta(nw),va(nw),vb(nw),s1(nw),s2(nw),s3(nw)
-
-k=(Cab*k_Cab+Car*k_Car+Cbrown*k_Brown+Cw*k_Cw+Cm*k_Cm)/N
-
-! ********************************************************************************
-! reflectance and transmittance of one layer
-! ********************************************************************************
-! Allen W.A., Gausman H.W., Richardson A.J., Thomas J.R. (1969), Interaction of
-! isotropic ligth with a compact plant leaf, Journal of the Optical Society of
-! American, 59:1376-1379.
-! ********************************************************************************
-
-! exponential integral: S13AAF routine from the NAG library
-
-where (k.le.0.0)
- tau=1
-end where
-where (k.gt.0.0.and.k.le.4.0)
- xx=0.5*k-1.0
- yy=(((((((((((((((-3.60311230482612224d-13 &
- *xx+3.46348526554087424d-12)*xx-2.99627399604128973d-11) &
- *xx+2.57747807106988589d-10)*xx-2.09330568435488303d-9) &
- *xx+1.59501329936987818d-8)*xx-1.13717900285428895d-7) &
- *xx+7.55292885309152956d-7)*xx-4.64980751480619431d-6) &
- *xx+2.63830365675408129d-5)*xx-1.37089870978830576d-4) &
- *xx+6.47686503728103400d-4)*xx-2.76060141343627983d-3) &
- *xx+1.05306034687449505d-2)*xx-3.57191348753631956d-2) &
- *xx+1.07774527938978692d-1)*xx-2.96997075145080963d-1
- yy=(yy*xx+8.64664716763387311d-1)*xx+7.42047691268006429d-1
- yy=yy-log(k)
- tau=(1.0-k)*dexp(-k)+k**2*yy
-end where
-where (k.gt.4.0.and.k.le.85.0)
- xx=14.5/(k+3.25)-1.0
- yy=(((((((((((((((-1.62806570868460749d-12 &
- *xx-8.95400579318284288d-13)*xx-4.08352702838151578d-12) &
- *xx-1.45132988248537498d-11)*xx-8.35086918940757852d-11) &
- *xx-2.13638678953766289d-10)*xx-1.10302431467069770d-9) &
- *xx-3.67128915633455484d-9)*xx-1.66980544304104726d-8) &
- *xx-6.11774386401295125d-8)*xx-2.70306163610271497d-7) &
- *xx-1.05565006992891261d-6)*xx-4.72090467203711484d-6) &
- *xx-1.95076375089955937d-5)*xx-9.16450482931221453d-5) &
- *xx-4.05892130452128677d-4)*xx-2.14213055000334718d-3
- yy=((yy*xx-1.06374875116569657d-2)*xx-8.50699154984571871d-2)*xx+9.23755307807784058d-1
- yy=exp(-k)*yy/k
- tau=(1.0-k)*dexp(-k)+k**2*yy
-end where
-where (k.gt.85.0)
- tau=0
-end where
-
-! transmissivity of the layer
-
-theta1=90.
-call tav_abs(theta1,refractive,t1)
-theta2=40.
-call tav_abs(theta2,refractive,t2)
-x1=1-t1
-x2=t1**2*tau**2*(refractive**2-t1)
-x3=t1**2*tau*refractive**2
-x4=refractive**4-tau**2*(refractive**2-t1)**2
-x5=t2/t1
-x6=x5*(t1-1)+1-t2
-r=x1+x2/x4
-t=x3/x4
-ra=x5*r+x6
-ta=x5*t
-
-! ********************************************************************************
-! reflectance and transmittance of N layers
-! ********************************************************************************
-! Stokes G.G. (1862), On the intensity of the light reflected from or transmitted
-! through a pile of plates, Proceedings of the Royal Society of London, 11:545-556.
-! ********************************************************************************
-
-delta=(t**2-r**2-1)**2-4*r**2
-beta=(1+r**2-t**2-sqrt(delta))/(2*r)
-va=(1+r**2-t**2+sqrt(delta))/(2*r)
-vb=sqrt(beta*(va-r)/(va*(beta-r)))
-s1=ra*(va*vb**(N-1)-va**(-1)*vb**(-(N-1)))+(ta*t-ra*r)*(vb**(N-1)-vb**(-(N-1)))
-s2=ta*(va-va**(-1))
-s3=va*vb**(N-1)-va**(-1)*vb**(-(N-1))-r*(vb**(N-1)-vb**(-(N-1)))
-RT(:,1)=s1/s3
-RT(:,2)=s2/s3
-
-end subroutine
diff --git a/src/models/tav_abs.f90 b/src/models/tav_abs.f90
deleted file mode 100644
index 9c03dea..0000000
--- a/src/models/tav_abs.f90
+++ /dev/null
@@ -1,60 +0,0 @@
-! ********************************************************************************
-! tav_abs.f90
-! ********************************************************************************
-! computation of the average transmittivity at the leaf surface within a given
-! solid angle. teta is the incidence solid angle (in radian). The average angle
-! that works in most cases is 40deg*pi/180. ref is the refaction index.
-! ********************************************************************************
-! Stern F. (1964), Transmission of isotropic radiation across an interface between
-! two dielectrics, Applied Optics, 3:111-113.
-! Allen W.A. (1973), Transmission of isotropic light across a dielectric surface in
-! two and three dimensions, Journal of the Optical Society of America, 63:664-666.
-! ********************************************************************************
-! version 5.02 (25 July 2011)
-! ********************************************************************************
-
-subroutine tav_abs(theta,refr,res)
-
-use dataSpec_P5B
-implicit none
-
-real(8), intent(in) :: theta, refr(nw)
-real(8), intent(out) :: res(nw)
-
-real(8) pi, thetarad
-real(8) refr2(nw)
-real(8) ax(nw),bx(nw),b0(nw),b1(nw),b2(nw)
-real(8) ts(nw),tp(nw),tp1(nw),tp2(nw),tp3(nw),tp4(nw),tp5(nw)
-
-pi=atan(1.)*4.
-thetarad=pi*theta/180.
-
-if (theta.eq.0.) then
- res=4.*refr/(refr+1.)**2
- return
-endif
-
-refr2=refr*refr
-ax=(refr+1.)**2/2.
-bx=-(refr2-1.)**2/4.
-
-if (thetarad.eq.pi/2.) then
- b1=0.
-else
- b1=dsqrt((sin(thetarad)**2-(refr2+1.)/2.)**2+bx)
-endif
-b2=sin(thetarad)**2-(refr2+1.)/2.
-b0=b1-b2
-ts=(bx**2/(6.*b0**3)+bx/b0-b0/2.)-(bx**2/(6.*ax**3)+bx/ax-ax/2.)
-tp1=-2.*refr2*(b0-ax)/(refr2+1.)**2
-tp2=-2.*refr2*(refr2+1.)*dlog(b0/ax)/(refr2-1.)**2
-tp3=refr2*(1./b0-1./ax)/2.
-tp4=16.*refr2**2*(refr2**2+1.)*dlog((2.*(refr2+1.)*b0-(refr2-1.)**2)/ &
- (2.*(refr2+1.)*ax-(refr2-1.)**2))/((refr2+1.)**3*(refr2-1.)**2)
-tp5=16.*refr2**3*(1./(2.*(refr2+1.)*b0-((refr2-1.)**2))-1./(2.*(refr2+1.) &
- *ax-(refr2-1.)**2))/(refr2+1.)**3
-tp=tp1+tp2+tp3+tp4+tp5
-res=(ts+tp)/(2.*sin(thetarad)**2)
-
-return
-end
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-grass/pktools.git
More information about the Pkg-grass-devel
mailing list