#DEFINITIONS:  -*-sh-*-
#
# Demonstration and test the use of CIA contiua in RT calculations, 
# together with other absorption tags.
#
# Based on Patrick Eriksson's TestClearSky2.arts
#
# Author: Stefan Buehler
# Date:   2013-03-08


Arts2 {

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

# (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 )

# Planck as blackbody radiation
Copy( blackbody_radiation_agenda, blackbody_radiation_agenda__Planck )

# on-the-fly absorption
Copy( propmat_clearsky_agenda, propmat_clearsky_agenda__OnTheFly )

# Include CIA in spectroscopy
Copy( abs_xsec_agenda, abs_xsec_agenda__withCIA )

# sensor-only path
Copy( ppath_agenda, ppath_agenda__FollowSensorLosPath )

# no refraction
Copy( ppath_step_agenda, ppath_step_agenda__GeometricPath )


# Number of Stokes components to be computed
#
IndexSet( stokes_dim, 1 )

# Clearsky = No scattering
#
cloudboxOff

# No jacobian calculation
#
jacobianOff


# Definition of species
# ---
# ATTENTION: You also have to include CIA in abs_xsec_agenda,
# otherwise you will not get any absorption for it!
abs_speciesSet( species=
            ["H2O-SelfContCKDMT100, H2O-ForeignContCKDMT100, H2O", 
             "O2",
             "N2-CIA-N2-0, N2-CIA-N2-1"] )

# Read CIA data
# ---
ReadXML(abs_cia_data, "spectroscopy/cia/hitran2011/hitran_cia2012_adapted.xml.gz")


# Line file and a matching coarse frequency grid
# ---
abs_linesReadFromSplitArtscat( 	
    abs_lines, 
    abs_species, 
    "spectroscopy/Perrin/", 
    1e9, 3000e9 )
VectorNLinSpace( f_grid, 100, 1e9, 3000e9 )


# A pressure grid rougly matching 0 to 80 km, in steps of 2 km.
# ---
VectorNLogSpace( p_grid, 41, 1000e2, 1 )


# Sort the line file according to species
# ---
abs_lines_per_speciesCreateFromLines


# Atmospheric scenario
# ---
AtmRawRead( t_field_raw, z_field_raw, vmr_field_raw, abs_species, 
            "testdata/tropical" )


# Weakly reflecting surface
# ---
VectorSetConstant( surface_scalar_reflectivity, 1, 0.8 )
Copy( surface_rtprop_agenda,
      surface_rtprop_agenda__Specular_NoPol_ReflFix_SurfTFromt_surface )


# Tb unit
# ---
StringSet( iy_unit, "PlanckBT" )


# Atmosphere and surface
# ---
AtmosphereSet1D
AtmFieldsCalc
Extract( z_surface, z_field, 0 )
Extract( t_surface, t_field, 0 )


# Definition of sensor position and LOS
# ---
VectorSet( rte_pos, [600e3] )
VectorSet( rte_los, [180] )
VectorSet( rte_pos2, [] )     # A dummy value

# Perform RT calculations
# ---
abs_xsec_agenda_checkedCalc
propmat_clearsky_agenda_checkedCalc
atmfields_checkedCalc
atmgeom_checkedCalc
cloudbox_checkedCalc
iyCalc

#WriteXML( "ascii", iy )
#WriteXML( "ascii", f_grid )
#WriteXML( "ascii", abs_species )


# Compare results to reference calculation:
MatrixCreate( iy_reference )
ReadXML( iy_reference,
         "TestRTwithCIA.iy_reference.xml")

Compare( iy, 
         iy_reference, 
         1e-10, 
	 "Brightness temperature with CIA N2 continuum differs from reference calculation" )

}