#
# A sample control file to demonstrate some of the RTE, WF and batch
# functions. You will need to edit some paths to get things working.
#
# Run this control file and type in AMI:
#    f=read_artsvar('full','f_mono');   
#    za=read_artsvar('full','za_pencil'); 
#    y=read_artsvar('full','y');   
#    plot(f/1e9,reshape(invrayjean(501e9,y),length(f),length(za))) 
#
# You should then get something similar to Figure 4 in Paper A of my thesis.
#
# By typing 
#    Y=read_artsvar('full','ybatch');  
# you get the result of the batch calculations.
#

# PE 2001-04-23
# The RTE and WF parts are updated while the batch part is not up to date.


### Spectroscopy: ############################################################
tgsDefine           { [
       "O3",
       "ClO",
       "H2O,H2O-ContStandardForeign,H2O-ContStandardSelf",
       "N2O",
       "N2-ContRosenkranzSelf" ] }
wfs_tgsDefine       { ["O3","ClO"] }
lineshapeDefine     { shape               = "Voigt_Kuntz6"
                      normalizationfactor = "quadratic"
                      cutoff              = -1 } 
#
linesReadFromHitran { filename = "/smiles_local/Qpack/Samples/Spectroscopy/lines501.4"
		      fmin     = 498e9
                      fmax     = 505e9 }
lines_per_tgCreateFromLines {}  
#
cont_descriptionInit{
}
cont_descriptionAppend{
    name       = "H2O-ContStandardSelf"
    parameters = [ 17.96e-34, 4.5 ]
}
cont_descriptionAppend{
    name       = "H2O-ContStandardForeign"
    parameters = [ 5.43e-35, 0.0 ]
}
cont_descriptionAppend{
    name       = "N2-ContRosenkranzSelf"
    parameters = [ 1.05e-38, 3.55 ]
}



### Load an atmosphere ########################################################
MatrixReadAscii (raw_ptz_1d) { "/smiles_local/arts/data/atmosphere/fascod/midlatitude-winter.tz.am" }
raw_vmrs_1dReadFromScenario  { "/smiles_local/arts/data/atmosphere/fascod/midlatitude-winter" }


### Monochromatic frequency grid ##############################################
VectorLinSpace (f_mono)    { start = 501.18e+9
                             stop  = 501.58e+9
		             step  = 1e6       }
VectorWriteAscii (f_mono) {""}


### The geoid #################################################################
# Set r_geoid to EARTH_RADIUS
r_geoidStd{}


### Vertical profiles #########################################################
VectorNLogSpace (p_abs)    { start = 1.013300e+05
                             stop  = 0.0914023919463313
                             n     = 101 }
AtmFromRaw1D               {}
h2o_absSet{}
n2_absSet{}
#hseSet{
#   pref  = 101330
#   zref  = 0
#   g0    = 9.81
#   niter = 2 }
hseFromBottom{
   g0    = 9.81
   niter = 2 }
#hseOff{}
hseCalc{}
VectorWriteAscii (p_abs)  {""}
VectorWriteAscii (t_abs)  {""}
VectorWriteAscii (z_abs)  {""}
#ArrayOfVectorWriteAscii (vmrs) {""}


### Calculate absorption ######################################################
absCalc           {}
abs_per_tgReduce  {}
#


### The ground ################################################################
#NumericSet (z_ground) { 0.0 }
#NumericSet (t_ground) { 290.0 }
#VectorSetLengthFromVector(e_ground,f_mono){value=0.75}
#groundOff{}
groundSet{
   z = 1e3
   e = 0.9
}


### Observation geometry ######################################################
NumericSet (z_plat) { 620e3 }
VectorNLinSpace (za_pencil) { start = 113
                              stop  = 114.1
                              n     = 11 }
#VectorSet  (za_pencil) { length = 1
#                         value  = 0 }
VectorWriteAscii( za_pencil ){""}
NumericSet (l_step) { 20e3 }


### Cosmic radiation ##########################################################
y_spaceStd{ "cbgr" }
#y_spacePlanck{300}


### Refraction ################################################################
#refrOff{}
refrSet{
  on    = 1
  model = "Boudouris"
  lfac  = 4
}
refrCalc{}


### Emission ##################################################################
emissionOn {}


### Line of sight (LOS) #######################################################
losCalc{}
VectorWriteBinary( z_tan ) {""}
sourceCalc{}
transCalc{}


### Calculate spectra #########################################################
yCalc{}
VectorWriteBinary( y ) {""}
#yTau{}




### Calculate Kx ##############################################################

### Init ###
#kxInit{}
#kbInit{}
kxAllocate( y ){
   ni = 4 
   nx = 224 }

### Species ###
absloswfsCalc{}
VectorCopy( k_grid, p_abs ){}
#kSpecies{    tag   = "O3"
#             unit  = "frac" }
kSpeciesAll{ unit = "frac" }
#kxAppend{}
kxPutInK{}

### Pointing off-set ###
kPointingOffSet{ delta = 0.001 }
#kxAppend{}
kxPutInK{}

### Temperature ###
VectorNLogSpace (k_grid)   { start = 1.013300e+05
                             stop  = 0.0914023919463313
                             n     = 21 }
#kTemp{}
kTempFast{}
#kxAppend{}
kxPutInK{}

### Save ###
MatrixWriteBinary(kx){""}  
ArrayOfIndexWriteBinary(kx_lengths){""}  
MatrixWriteBinary(kx_aux){""}  
ArrayOfStringWriteBinary(kx_names){""}

Exit{}


### Conversion to different temperature scales ################################
VectorCopy(f_sensor,f_mono){}
VectorCopy(za_sensor,za_pencil){}
yTRJ{}
#yTB{}
#VectorPlanck(y_cal1,f_sensor){temp=78}
#VectorPlanck(y_cal2,f_sensor){temp=300}
#VectorSetLengthFromVector(i_cal1,y_cal1){value=78}
#VectorSetLengthFromVector(i_cal2,y_cal2){value=300}          
#yLoadCalibration{}


### Save spectra #############################################################
VectorWriteBinary( y ) {""}



### Produce batch data for varoius variables ##################################

StringSet( batchname ) { "fullbatch" }

sDiagonalLengthFromVector ( p_abs )
{ 
  stddev     = 3
}

BatchdataGaussianTemperatureProfiles{
   n     = 20
}

sDiagonalLengthFromVector ( p_abs )
{ 
  stddev     = 0.3
}

BatchdataGaussianSpeciesProfiles{
   n       = 20
   do_tags = ["O3-666"]
   unit    = "frac"
}


BatchdataUniformOffSets ( za_sensor ) {
   n    = 20
   low  = -0.01
   high = 0.01
}

# These lines produce thermal noise, but the covariance cannot be 
# created for the selected f_sensor and za_sensor without 
# diagonal/banded covariance matrices.
#
#sDiagonalLengthFromVectors ( f_sensor, za_sensor )
#{ 
#  stddev     = 0.2
#}
#
#BatchdataGaussianZeroMean {
#   n       = 20
#   varname = "noise"
#}

BatchdataSinusoidalRippleNoCorrelations( f_sensor, za_sensor )
{
  n         = 20
  period    = 20e6
  amplitude = 0.2
  pdf       = "gaussian"
  varname   = "ripple"
}


### Calculate spectra for the batch data,  ####################################
### convert to RJ and add noise and ripple

ybatchAbsAndRte{
   ncalc     = 10
   do_t      = 1
   t_file    = ""
   do_z      = 1
   z_file    = ""
   do_tags   = ["O3-666"]
   tag_files = [""]
   do_f      = 0
   f_file    = ""
   do_za     = 1
   za_file   = ""
}

MatrixTRJ( ybatch, ybatch ){}
#MatrixFillWithVector( ybatch, y ){ 10 }

ybatchAdd{ "ripple" }
#ybatchAdd{ "noise" }

MatrixWriteBinary( ybatch ) {""}