Thank you Freddy! This makes sense.
I tried with flat reflective surface and lambertian surface with
absorption only. I extended the 1D atmosphere to 3D so I can provide
surface temperature (set it to a really small value) but in theory
since there is no scattering the path that contributes is unchanged. I
found that the flat reflective surface result is correct but the
lambertian one is much smaller than the expected pi ratio. I tried to
print ppath and I found los changing quite a lot. I am wondering if
this is the expected behavior and maybe there is something about the
3D setup I did not understand properly. I copied the main part of the
code below.
Thanks!
Mengqi
@arts_agenda
def propmat_clearsky_agenda(ws):
ws.Ignore(ws.rtp_mag)
ws.Ignore(ws.rtp_los)
ws.propmat_clearskyInit()
ws.propmat_clearskyAddConts()
ws.propmat_clearskyAddLines()
@arts_agenda
def gas_scattering_agenda(ws):
ws.Ignore(ws.rtp_vmr)
ws.gas_scattering_coefXsecConst(ConstXsec=4.65e-31)
ws.gas_scattering_matIsotropic()
# surface scattering agenda
# lambertian
@arts_agenda
def iy_surface_agenda(ws):
ws.iySurfaceInit()
ws.Ignore(ws.dsurface_rmatrix_dx)
ws.Ignore(ws.dsurface_emission_dx)
ws.iySurfaceLambertian()
ws.iySurfaceLambertianDirect()
# flat reflective surface
# @arts_agenda
# def iy_surface_agenda(ws):
# ws.iySurfaceInit()
# ws.iySurfaceFlatReflectivity()
# ws.iySurfaceFlatReflectivityDirect()
# generate atmosphere data
dataset_path =
'/home/mandy/Github/arts/build_new/afgl_1986-us_standard.nc'
save_path = '/home/mandy/Github/arts/controlfiles/testdata/'
data = generate_atmos_arts(dataset_path, save_path)
#
=============================================================================
# open workspace
#
=============================================================================
ws = Workspace()
ws.verbositySetScreen(level=2)
#
=============================================================================
# generate atmosphere data
#
=============================================================================
dataset_path =
'/home/mandy/Github/arts/build_new/afgl_1986-us_standard.nc'
save_path = '/home/mandy/Github/arts/controlfiles/testdata/'
data = generate_atmos_arts(dataset_path, save_path)
ws.ReadHITRAN(filename='/home/mandy/Github/MiAtmosphere/HITRAN/ALL.par',
hitran_type="Online", abs_lines=ws.abs_lines)
#
=============================================================================
# select/define agendas
#
=============================================================================
ws.LegacyContinuaInit()
ws.PlanetSet(option="Earth")
# cosmic background radiation
ws.iy_space_agendaSet( option="CosmicBackground" )
# sensor-only path
ws.ppath_agendaSet( option="FollowSensorLosPath" )
# no refraction
ws.ppath_step_agendaSet( option="GeometricPath" )
# main agenda
ws.iy_main_agendaSet( option="Clearsky")
# water agenda
ws.water_p_eq_agendaSet()
# surface agenda
ws.iy_surface_agenda = iy_surface_agenda
ws.ArrayOfStringSet( ws.iy_aux_vars,
[ "Optical depth",
"Radiative background"
] )
ws.propmat_clearsky_agenda=propmat_clearsky_agenda
# gas scattering agenda
ws.gas_scattering_agenda = gas_scattering_agenda
ws.NumericSet( ws.ppath_lmax, 1e10)
#
=============================================================================
# basic conditions
#
=============================================================================
# Postion and line-of-sight of sensor
sensor_pos = np.array([[600e+3, 0., 0.]])
# Sensor looking direction in zenith angle (0 = upwards, 180 =
downward) and
# azimuth angle ( 0 = North, 90 = east)
sensor_los = np.array([[180, 0]])
ws.sensor_pos = sensor_pos
ws.sensor_los = sensor_los
ws.VectorSet(ws.rte_pos2, [])
# define environment
#
=============================================================================
# Number of Stokes components to be computed
ws.IndexSet(ws.stokes_dim, 1)
# Read the spectroscopic line data from the ARTS catalogue and
# create the workspace variable `lines'.
ws.ReadHITRAN(filename='/home/mandy/Github/MiAtmosphere/HITRAN/ALL.par',
hitran_type="Online", abs_lines=ws.abs_lines)
ws.abs_linesNormalization(ws.abs_lines, "VVH")
# Frequency grid
c = 299792458
# Set frequency
#wavelengths = np.linspace(1250e-9, 1177e-9, 74)
#wavelengths = np.linspace(700e-9, 400e-9, 74)
wavelengths = np.linspace(1178e-9, 1177e-9, 2)
f_grid = c / wavelengths
ws.f_grid = f_grid
####
# set a simple blackbody sun
ws.Touch(ws.stars)
ws.starsAddSingleBlackbody(distance=1.495978707e11, latitude=0.,
longitude=0.)
ws.Print(ws.stars, 2)
# Reference ellipsoid
ws.refellipsoidEarth(ws.refellipsoid, "Sphere")
# A pressure grid rougly matching 0 to 80 km, in steps of 2 km.
#ws.p_grid = data.p.to_numpy()
ws.p_grid = np.array([1.01040e+05, 2.54000e-03]) # 0 and 120km
# Atmospheric dimensionality and lat/lon grids
nlat = 3
nlon = 5
ws.VectorNLinSpace(ws.lat_grid, nlat, -90., 90.)
ws.VectorNLinSpace(ws.lon_grid, nlon, -180., 180.)
ws.AtmosphereSet3D()
# Definition of species
ws.abs_speciesSet(species=
["H2O", "O3", "N2O", "CO", "CH4", "CO2", "O2"])
# This separates the lines into the different tag groups and creates
# the workspace variable `abs_lines_per_species':
ws.abs_lines_per_speciesCreateFromLines()
# Load atmospheric data
ws.AtmRawRead(basename="../controlfiles/testdata/afglUS")
ws.propmat_clearsky_agendaAuto()
ws.AtmFieldsCalcExpand1D()
#Get ground altitude (z_surface) from z_field
ws.MatrixSetConstant(ws.z_surface, nlat, nlon, 0.)
ws.ArrayOfStringSet(ws.surface_props_names, ["Skin temperature"])
#ws.Tensor3SetConstant(ws.surface_props_data, 1, nlat, nlon,
ws.t_field.value[0,0,0])
ws.Tensor3SetConstant(ws.surface_props_data, 1, nlat, nlon, 1e-16)
# #print(ws.t_field.value[0,0,0])
# Set surface relectivity
#ws.surface_reflectivity = np.array([[[1.]]])
ws.surface_scalar_reflectivity = [1]
# No jacobian calculations
ws.jacobianOff()
# No particulate scattering
ws.cloudboxOff()
# No sensor model
ws.sensorOff()
ws.StringSet( ws.iy_unit, "1" )
ws.Extract( ws.z_surface, ws.z_field, 0 ) # z_field is all the
altitude and z_surface is altitude at the surface
ws.Extract( ws.t_surface, ws.t_field, 0 ) # t_field is all the temperature
ws.WriteXML( ws.output_file_format, ws.z_field )
ws.WriteXML( ws.output_file_format, ws.t_field )
#Switch off gas scattering
ws.IndexSet(ws.gas_scattering_do, 0)
#Switch on stars
ws.IndexSet(ws.stars_do, 1)
# Check model atmosphere
ws.propmat_clearsky_agendaAuto()
ws.propmat_clearsky_agenda_checkedCalc()
ws.atmfields_checkedCalc()
ws.atmgeom_checkedCalc()
ws.cloudbox_checkedCalc()
ws.sensor_checkedCalc()
ws.lbl_checkedCalc()
# the actual simulation
ws.yCalc()
# output iy file
ws.WriteXML(ws.output_file_format, ws.y)
------------------------------------------------------------------------
*From:* Manfred Brath <manfred.br...@uni-hamburg.de>
*Sent:* Monday, October 17, 2022 5:56:40 PM
*To:* Xia Mengqi; arts_users.mi@lists.uni-hamburg.de
*Subject:* Re: Fwd: [arts-users] Direct radiation is subtracted when
it should not
Hello Mengqi,
Am 17.10.22 um 16:33 schrieb Xia Mengqi:
Hi Freddy,
Thank you so much and this is very helpful! Just to make sure I fully
understand how to use it in the code -- According to ARTS output, it
seems I need to provide the skin temperature as well. I noticed that
in one provided example it has "
ws.ArrayOfStringSet(ws.surface_props_names, ["Skin temperature"])
ws.Tensor3SetConstant(ws.surface_props_data, 1, nlat, nlon,
ws.t_field.value[0, 0, 0])"
I'm wondering if there is a simple way to do this for 1D.
so far you cannot use a 1d atmosphere for simulation with a direct
source like the sun because the atmosphere dimension and (radiative
transfer) geometry are coupled in ARTS. Since you have a zenith and an
azimuth dependency for simulation with a direct source, you have to
use a 3d atmosphere in ARTS as simulation of 1d atmospheres have only
zenith dependency in ARTS.
Also, is it possible to just not include surface emission?
No, but you can set the surface temperature to a small value greter
than 0 K. This would have a similar effect.
Cheers,
Freddy
--
----------------------------------------------------------
Dr. Manfred Brath
Radiation and Remote Sensing
Meteorological Institute
Universität Hamburg
Bundesstraße 55
D-20146 Hamburg
Room 1535
Tel: +49 40 42838-8786
