iy_01.arts

Arts2 {

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

# Basic settings of simulation
    AtmosphereSet1D
    IndexSet( stokes_dim, 1 )
    StringSet( iy_unit, "PlanckBT" )

# monochromatic frequency grid
    #VectorLinSpace( f_grid, 22e9, 70e9, 1e7 )
    #VectorNLinSpace( f_grid, 500, 10e9, 230e9 )
    VectorNLinSpace( f_grid, 1125, 5e9, 500e9 )



# Common microwave sensor settings
    output_file_formatSetZippedAscii

    # 1. General
    # modify the maximum propagation step, from the default to 250 m :
    NumericSet( ppath_lmax, 250 )

    # Agenda for scalar gas absorption calculation
    AgendaSet( abs_xsec_agenda ){
      Ignore(abs_nlte)
      abs_xsec_per_speciesInit
      abs_xsec_per_speciesAddPredefinedO2MPM2020
      abs_xsec_per_speciesAddConts
    }
    
    # Surface
    Copy( surface_rtprop_agenda,
    surface_rtprop_agenda__Specular_NoPol_ReflFix_SurfTFromt_surface )

    # (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 propmat_clearsky_agenda to use on-the-fly absorption
    #Copy( propmat_clearsky_agenda, propmat_clearsky_agenda__OnTheFly )
    AgendaSet( propmat_clearsky_agenda ){
      Ignore(rtp_mag)
      Ignore(rtp_los)
      propmat_clearskyInit
      propmat_clearskyAddXsecAgenda
      propmat_clearskyAddLines
      propmat_clearskyForceNegativeToZero
    }
    
    # Clearsky = No scattering
    #
    cloudboxOff

    # No jacobian calculation
    #
    jacobianOff

    # No sensor properties
    # ---
    sensorOff
    
# ====================================================================

# Spectroscopy
    abs_speciesSet( species=[ "H2O, H2O-SelfContCKDMT252, H2O-ForeignContCKDMT252",
                              "O2-MPM2020",                              
                              "N2,  N2-CIAfunCKDMT252, N2-CIArotCKDMT252",
                              "O3" ] )

    ReadARTSCAT( abs_lines, "instruments/metmm/abs_lines_metmm.xml.gz" )
    abs_linesSetCutoff(option="ByLine", value=750e9)
    abs_linesSetNormalization(option="VVH")
    abs_lines_per_speciesCreateFromLines
    abs_lines_per_speciesSetCutoffForSpecies(option="ByLine", value=5e9, species_tag="O3")
#    WriteXML( "zascii", abs_lines, "instruments/metmm/abs_lines_metmm.xml.gz" )


# ====================================================================

# Set surface reflectivity
# Reflectivity = 0.4; emissivity = 0.6
    VectorSetConstant( surface_scalar_reflectivity, 1, 0.05 )

# Atmospheric profiles
    ReadXML( batch_atm_fields_compact, "testdata/garand_profiles.xml.gz" )

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

# Absorption lookup table
    abs_xsec_agenda_checkedCalc
    lbl_checkedCalc

# Setting the agenda for batch calculation 
# Garand profiles have 42 different. We will make RT calculations for all of them.
    ArrayOfMatrixCreate( out )
    #Tensor4Create( temp_tensor )
    ArrayOfMatrixCreate( temp_tensor )
    VectorCreate( temp_vector )
    MatrixCreate( temp_matrix )
    
    AgendaSet( ybatch_calc_agenda ){
  
# Extract the atmospheric profiles for this case:
    Extract(
    atm_fields_compact, 
    batch_atm_fields_compact, 
    ybatch_index
    )

# Split up *atm_fields_compact* to
# generate p_grid, t_field, z_field, vmr_field:
    AtmFieldsAndParticleBulkPropFieldFromCompact

# Optionally set Jacobian parameters.
    jacobianOff

# No scattering
    cloudboxOff
    
# No sensor
    sensorOff
    
# Definition of sensor position and LOS
# ---
    MatrixSetConstant( sensor_pos, 1, 1, 850e3 )
    MatrixSet( sensor_los, [ 180 ] )
    VectorSet( rte_pos, [850e3] )
    VectorSet( rte_los, [180] )
    VectorSet( rte_pos2, [] )   

# get some surface properties from corresponding atmospheric fields
    Extract( z_surface, z_field, 0 )
    Extract( t_surface, t_field, 0 )

# Checks
    atmfields_checkedCalc
    atmgeom_checkedCalc
    cloudbox_checkedCalc
    sensor_checkedCalc
# Perform RT calculations
    ArrayOfStringSet( iy_aux_vars,
     [
      "Absorption, summed", 
      "Absorption, species 0",
      "Absorption, species 1",
      "Absorption, species 2",
      "Absorption, species 3",
      "iy", 
      "Transmission",
      "Optical depth"
     ] )

    iyCalc
    
    #=====================================================================
            #### Output ####
    #=====================================================================
      
    # ybatchCalc braucht y, y_aux und jacobian als Output. Da diese nicht
    # erzeugt werden, kann ARTS mit Touch() mitgeteilt werden, dass die 
    # Variable benutzt wurde.
    Touch( y )
    Touch( y_aux )
    Touch( jacobian )

    # Zusammenfassen aller Output-Variablen in die Variable out

    Matrix1ColFromVector( temp_matrix, f_grid )
    Append( out, temp_matrix )

    Matrix1ColFromVector( temp_matrix, p_grid )
    Append( out, temp_matrix )

    Reduce( temp_vector, z_field )
    Matrix1ColFromVector( temp_matrix, temp_vector )
    Append( out, temp_matrix )

    Reduce( temp_vector, t_field )
    Matrix1ColFromVector( temp_matrix, temp_vector )
    Append( out, temp_matrix )

    Print(iy_aux)
    Exit
    Extract( temp_tensor, iy_aux, 0 ) 
    Reduce( temp_vector, temp_tensor )
    Matrix1ColFromVector( temp_matrix, temp_vector )
    Append( out, temp_matrix )

    Extract( temp_tensor, iy_aux, 1 ) 
    Reduce( temp_matrix, temp_tensor )
    Append( out, temp_matrix )

    Extract( temp_tensor, iy_aux, 2 ) 
    Reduce( temp_matrix, temp_tensor )
    Append( out, temp_matrix )

    Extract( temp_tensor, iy_aux, 3 )
    Reduce( temp_matrix, temp_tensor )
    Append( out, temp_matrix )

    Extract( temp_tensor, iy_aux, 4 ) 
    Reduce( temp_vector, temp_tensor )
    Matrix1ColFromVector( temp_matrix, temp_vector )
    Append( out, temp_matrix )

    Extract( temp_tensor, iy_aux, 5 ) 
    Reduce( temp_matrix, temp_tensor )
    Append( out, temp_matrix )

    Extract( temp_tensor, iy_aux, 6 ) 
    Reduce( temp_matrix, temp_tensor )
    Append( out, temp_matrix )

    Extract( temp_tensor, iy_aux, 7 )
    Reduce( temp_matrix, temp_tensor )
    Append( out, temp_matrix )
    
    Append( out, iy )

    # out = ( f_grid, p_grid, z_field, t_field, vmr H20, 
    # abs H20, abs O2, abs[f,z], PlacnkBT )
    WriteNetCDFIndexed( ybatch_index, out, "Output/2021-07-18_2322_out" )
}

# How many simulations do we want to perform?
# All atmospheres, or manually set the number (uncomment IndexSet line)
    nelemGet( ybatch_n, batch_atm_fields_compact )
    #IndexSet(ybatch_n, 1)
# ====================================================================

# Execute the batch calculations:
# First check, then execute the batch RT calculations
    propmat_clearsky_agenda_checkedCalc
    ybatchCalc( robust=1 )
# ====================================================================

# Store results
    WriteXML( "ascii", f_grid, "Output/2021-07-18_2322_fgrid.xml" )
    WriteXML( "ascii", ybatch_n, "2021-07-18_2322_out_ybatch_n.xml" )
    
}