#DEFINITIONS:  -*-sh-*-
#
# Demonstration and test of a single radio link calculations, using iyCalc.
#
# The link between a satellite and a surface point is simulated. The refracted
# path is at about 5 degress above the horizon (at the surface point). 
#
# 2012-08-21, Patrick Eriksson

Arts2 {

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

# Agenda for scalar gas absorption calculation
Copy(abs_xsec_agenda, abs_xsec_agenda__noCIA)

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


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

# Frequency grid 
#
VectorSet( f_grid, [ 1e9, 5e9, 10e9 ] )


# A pressure grid rougly matching 0 to 80 km, in steps of 250.
#
VectorNLogSpace( p_grid, 321, 1013e2, 1 )


# Definition of species
# 
abs_speciesSet( species=
                ["H2O-PWR98", 
                 "N2-SelfContStandardType",
                 "O2-PWR93"] )


# No line data needed here
# 
abs_lines_per_speciesSetEmpty


# Dimensionality of the atmosphere
#
AtmosphereSet1D

# Atmospheric profiles
# 
AtmRawRead( t_field_raw, z_field_raw, vmr_field_raw, abs_species, 
            "testdata/tropical" )
#
AtmFieldsCalc( t_field, z_field, vmr_field, p_grid, 
               lat_grid, lon_grid, t_field_raw, z_field_raw, 
               vmr_field_raw, atmosphere_dim, 3 )


# Surface altitude
MatrixSetConstant( z_surface, 1, 1, 0 )


# No scattering
#
cloudboxOff

# No jacobian calculations
#
jacobianOff

# Check model atmosphere
# 
atmfields_checkedCalc
atmgeom_checkedCalc
cloudbox_checkedCalc


# Propagation path agendas and variables
#

# transmitter-receiver path
Copy( ppath_agenda, ppath_agenda__TransmitterReceiverPath )

# refracted path
Copy( ppath_step_agenda, ppath_step_agenda__RefractedPath )
Copy( refr_index_air_agenda, refr_index_air_agenda__GasThayer  )
#
NumericSet( ppath_lmax, 10e3 )
NumericSet( ppath_lraytrace, 100 )


# Radiative transfer agendas and variables
#
VectorSet( rte_los, [] )   # Dummy value
#
AgendaSet( iy_transmitter_agenda ){
  Ignore( rtp_pos )
  Ignore( rtp_los )
  MatrixUnitIntensity( iy, stokes_dim, f_grid )
}
# Here we use the numerical defocusing method
AgendaSet( iy_main_agenda ){
  Ignore( rte_los )
  iyRadioLink
}


# Postion of sensor/receiver and transmitter
#
VectorSet( rte_pos, [ 80e3 ] )
VectorSet( rte_pos2, [ 0, 5.1 ] )


# Auxilary variables
#
ArrayOfStringSet( iy_aux_vars, [ 
    "Pressure", "Temperature", "VMR, species 0", 
    "Atmospheric loss", "Absorption, summed", 
    "Free space loss", "Free space attenuation", 
    "Defocusing loss", "Extra path delay", "Bending angle" ] )


abs_xsec_agenda_checkedCalc
propmat_clearsky_agenda_checkedCalc

# Run and save
#
iyCalc


# Set "Atmospheric loss" is main variable:
# (If yCalc is used, iyReplaceFromAux shall be included in iy_main_agenda)
#
iyReplaceFromAux( iy, iy_aux, iy_aux_vars, jacobian_do, "Atmospheric loss" )


#WriteXML( "ascii", iy, "iyREFERENCE.xml" )
#WriteXML( "ascii", ppath, "ppath.xml" )


# OK?
# ---
MatrixCreate( iyREFERENCE )
ReadXML( iyREFERENCE, "iyREFERENCE.xml" )
Compare( iy, iyREFERENCE, 1e-5)

}