#DEFINITIONS:  -*-sh-*-

# This simple control file can be used to understand the effects of 
# scattering on the radiation. It uses DOIT to compute the 
# radiative calculations.
#
# Author: Jakob Doerr
# 

Arts2 {

INCLUDE "general/general.arts"
INCLUDE "general/continua.arts"
INCLUDE "general/agendas.arts"
INCLUDE "general/agendasDOIT.arts"
INCLUDE "general/planet_earth.arts"

# 1.General Settings:---------------------------------------------
#----------------------------------------------------------------- 

# Agenda for scalar gas absorption calculation
Copy(abs_xsec_agenda, abs_xsec_agenda__noCIA)
# (standard) emission calculation
Copy( iy_main_agenda, iy_main_agenda__Emission )
# cosmic background radiation
Copy( iy_space_agenda, iy_space_agenda__CosmicBackground )
# standard surface agenda (i.e., make use of surface_rtprop_agenda)
Copy( iy_surface_agenda, iy_surface_agenda__UseSurfaceRtprop )
# sensor-only path
Copy( ppath_agenda, ppath_agenda__FollowSensorLosPath )
# no refraction
Copy( ppath_step_agenda, ppath_step_agenda__GeometricPath )

# Set out put file format
#------------------------
output_file_formatSetAscii

# Define f_grid
#--------------
VectorSet( f_grid, [229.0e9] )

#Set stokes dim
#--------------
IndexSet (stokes_dim, 1)

#def of atmosphere
#-----------------
AtmosphereSet1D

#No jacobian calculations
#-----------------
jacobianOff

# No sensor
sensorOff

# Modifiy the maximum propagation step, from the default(10.e3)
# to 250 m:---------------------------------------------------
NumericSet( ppath_lmax, 250 )

# 2. Prepare atmospheric (clearsky) data------------------------------
#---------------------------------------------------------------------
# Set absorption species
abs_speciesSet( species=[ "H2O-PWR98",
                          "O3",
                          "O2-PWR93",
                          "N2-SelfContStandardType"
                        ] )

# Read atmospheric data
ReadXML(batch_atm_fields_compact,"testdata/chevallierl91_all_extract.xml")
Extract(atm_fields_compact,batch_atm_fields_compact,1)

# add constant profiles for O2 and N2
atm_fields_compactAddConstant( name="abs_species-O2", value=0.2095 )
atm_fields_compactAddConstant( name="abs_species-N2", value=0.7808 )

# 3. Absorption-------------------------------------------------
#---------------------------------------------------------------
AtmFieldsFromCompact
atmfields_checkedCalc(bad_partition_functions_ok=1)
# Read Catalog (needed for O3):
abs_linesReadFromSplitArtscat( basename="spectroscopy/Perrin_artscat34/",
                                fmin=0e9, fmax=1000e9 )

abs_lines_per_speciesCreateFromLines

abs_lookupSetup
abs_xsec_agenda_checkedCalc
abs_lookupCalc

# absorption from LUT
Copy( propmat_clearsky_agenda, propmat_clearsky_agenda__LookUpTable )
propmat_clearsky_agenda_checkedCalc


# 4. Surface properties-----------------------------------------
#---------------------------------------------------------------
# Set surface reflectivity (=1-emissivity)
# corresponds to emissivity=0.75
VectorSetConstant( surface_scalar_reflectivity, 1, 0. )

Copy( surface_rtprop_agenda,
      surface_rtprop_agenda__Specular_NoPol_ReflFix_SurfTFromt_surface )


# 5.Specification of cloud-----------------------------------------------
# -----------------------------------------------------------------------
ScatSpeciesInit
ScatElementsPndAndScatAdd(
  scat_data_files=["scattering/H2O_ice/MieSphere_R1.50000e+02um.xml"],
  					   pnd_field_files=["./input/pndfield_input.xml"])

scat_dataCalc
scat_data_checkedCalc										  

# 6. The calculations ------------------------------------------------------
# ----------------------------------------------------------------------
AtmFieldsFromCompact

cloudboxSetManually(p1=101300,p2=1000,lat1=0,lat2=0,lon1=0,lon2=0)

pnd_fieldCalcFrompnd_field_raw
# If you want more scattering elements, you can enlarge them here by a factor:
Tensor4Scale(pnd_field,pnd_field,4)


# Get case-specific surface properties from corresponding atmospheric fields
Extract( z_surface, z_field, 0 )
Extract( t_surface, t_field, 0 )

  # Consistency checks
atmfields_checkedCalc( bad_partition_functions_ok = 1 )
atmgeom_checkedCalc
cloudbox_checkedCalc  

AgendaSet( doit_conv_test_agenda ){
  doit_conv_flagAbsBT( epsilon=[0.01],
                       max_iterations=100,
                       nonconv_return_nan=1) 
  Print( doit_iteration_counter, 0 )

}

doit_za_interpSet( interp_method="linear" )
DOAngularGridsSet( doit_za_grid_size, scat_aa_grid, scat_za_grid,
                   40, 20, "" )

AgendaSet( doit_mono_agenda ){
  # Prepare scattering data for DOIT calculation (Optimized method):
  Ignore(f_grid)
  DoitScatteringDataPrepare
  # Perform iterations: 1. scattering integral. 2. RT calculations with 
  # fixed scattering integral field, 3. convergence test 
  doit_i_field_monoIterate(accelerated=1)
  }
   
DoitInit
DoitGetIncoming( rigorous=0 )
doit_i_fieldSetClearsky
DoitCalc

WriteXML("ascii",doit_i_field, "./results/DOIT_ifield.xml")
WriteXML("ascii",f_grid, "./results/f_grid.xml")
WriteXML("ascii",p_grid, "./results/p_grid.xml")
WriteXML("ascii",scat_za_grid, "./results/za_grid.xml")
WriteXML("ascii",cloudbox_limits, "./results/cloudbox_limits.xml")
WriteXML("ascii",t_field, "./results/t_field.xml")
WriteXML("ascii",vmr_field, "./results/vmr_field.xml")
WriteXML("ascii",pnd_field, "./results/pnd_field.xml")

Tensor4Scale(pnd_field,pnd_field,0)
DoitCalc
WriteXML("ascii",doit_i_field, "./results/DOIT_ifield_clearsky.xml")
} # End of Main