#DEFINITIONS:  -*-sh-*-
#
# Demonstration and test of a single radio link calculation, using yCalc.
#
# The link between two satellites is simulated.
#
# 2012-08-20, 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, [ 0.5e9, 1e9 ] )


# 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
#
AgendaSet( iy_transmitter_agenda ){
  Ignore( rtp_pos )
  Ignore( rtp_los )
  MatrixUnitIntensity( iy, stokes_dim, f_grid )
}
AgendaSet( iy_main_agenda ){
  Ignore( rte_los )
  iyRadioLink( defocus_method = 2, defocus_shift = 1e-3 ) 
}



# Sensor/receiver and transmitter
#
MatrixSet( sensor_pos, [ 600e3 ] )
MatrixSet( sensor_los, [     0 ] )    # Dummy value
MatrixSet( transmitter_pos, [ 600e3, 47.8 ] )
#
sensorOff
sensor_checkedCalc


# Auxilary variables
#
# These are the ones that can be obtained by yCalc.
#
ArrayOfStringSet( iy_aux_vars, [ 
  "Bending angle", 
  "Extra path delay",
  "Impact parameter" ] )



# Expected results
#
VectorCreate( yREFERENCE )
NumericCreate( bending_angleREFERENCE )
NumericCreate( extra_pdREFERENCE )
NumericCreate( impact_parREFERENCE )
#
# Values set 2013-10-07, on Patrick's desktop (Ubunt 13.04)
VectorSet( yREFERENCE, [1.07137e-15, 1.0652e-15 ] )
NumericSet( bending_angleREFERENCE, 0.377976 )
NumericSet( extra_pdREFERENCE, 259.05e-09 )
NumericSet( impact_parREFERENCE, 6.38905e+06 )
#
# Free space, defocusing and atmospheric loss part of y, and no need to check
# those terms 

abs_xsec_agenda_checkedCalc
propmat_clearsky_agenda_checkedCalc


# Calculate and check iy
#
yCalc
#
#Print(y)
#Print(y_aux)
#
Compare( y, yREFERENCE, 2e-17 )


# Bending angle
#
VectorCreate( vdummy )
NumericCreate( bending_angle )
Extract( vdummy, y_aux, 0 )
Extract( bending_angle, vdummy, 0 )
Compare( bending_angle, bending_angleREFERENCE, 1e-4 )


# Extra path delay
#
NumericCreate( extra_pd )
Extract( vdummy, y_aux, 1 )
Extract( extra_pd, vdummy, 0 )
Compare( extra_pd, extra_pdREFERENCE, 1e-9 )


# Impact parameter
#
NumericCreate( impact_par )
Extract( vdummy, y_aux, 2 )
Extract( impact_par, vdummy, 0 )
Compare( impact_par, impact_parREFERENCE, 5 )

}