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