ImaAtCor: A Physically Based Tool for Combined Atmospheric and Topographic Corrections of Remote Sensing Images

: ImaAtCor is a tool for the simultaneous correction of remotely sensed data from atmospheric and topographic effects, including second-order corrections, such as adjacency effects. The implemented approach is physically based and uses MODTRAN for computation of the main radiometric quantities. A user-friendly, comprehensive, and exhaustive graphic interface allows the user to choose from different correction levels. The various panels allow one to set all the parameters to correctly characterize the atmosphere and define the acquisition and illumination geometries. The tool provides a number of facilities to easily manage the correction process for a wide range of sensor data, including the ability to process multiple data in batch mode, which is very useful for dealing with temporal series. Under the inclusion of topographic correction, this tool allows the user to select a digital elevation model that is automatically resampled to the image resolution and processed to obtain the parameters for radiometric transfer modeling and the subsequent harmonization of the surface with the model inversion. This tool also includes utilities for the pre-processing of PRISMA data.


Introduction
The increased performance of satellite sensors requires a high-quality pre-processing chain to support remote sensing analysis in many fields [1]. In this context, atmospheric correction is a critical step to maximize information content [2][3][4][5]. In contexts of terrain roughness, remotely sensed data must also be corrected for topographic effects because the irradiation on a slope strongly varies with the relative angle between the sun and the slope azimuth, thus impacting the remotely sensed signal [6][7][8]. This effect has to be properly modeled and corrected to enhance image quality. As sensor performance has increased, empirical approaches for atmospheric and topographic correction have been gradually replaced with physically based procedures.
In the literature, physically based models for atmospheric correction are often matched with empirical ones for topographic correction in the case of terrain roughness [9,10]. Although these studies have shown good results, a number of works have demonstrated the advantages of using a unique physically based model for the simultaneous correction of both atmospheric and topographic effects [11][12][13][14][15].
In this paper, the ImaAtCor (Image Atmospheric Correction) software is presented with an exhaustive and user-friendly GUI (graphic user interface); this tool's facilities can easily manage the atmospheric correction process. This software allows one to atmospherically and topographically correct a wide range of remote sensing data and implements a fully physically based approach based on MODTRAN 4 simulations. The software, under the assumption of a Lambertian surface, can also implement a model for the simulation and correction of the impact of the environment on a remotely sensed signal. The models implemented in the tool were tested in various contexts and showed great potential for improving remotely sensed products, as shown in [4]. Many facilities to manage multiple or redundant applications are provided to reuse the calculated quantities and speed up the processing activities. After the launch of the PRISMA mission, the software was equipped with a set of features for the pre-processing and re-calibration of PRISMA data.

Atmospheric Model
This software implements different models with different complexity levels. When both adjacency and topographic correction are selected, the model is represented by the following equation from [4]: where L is the sensor radiance; Lpath is the path radiance (the solar signal reflected by the atmosphere and not reaching the surface); is the target reflectance; is the background reflectance; S is the spherical albedo; and A, Ad, and B are parameters that do not depend on the surface properties and can be expressed as a function of the downward and upward (diffuse and direct) transmittances, of the solar zenith angle and solar flux at the top of the atmosphere. In particular, A does not depend on upward diffuse transmittance, B does not depend on direct upward transmittance, and Ad only depends on direct transmittance [4]. The first term of the second member represents the contribution of the target due to the direct light that is directly reflected onto the sensor. The impact of the sloped terrain on this term is modelled by considering the angle , defined as the angle between the incoming solar radiation and the normal to the surface that, for an inclined surface, is different from the sun zenith angle [13]. For angles above π/2, is set to π/2. Under this condition, the target is in the shadow, and the term is zero. The second term accounts for the diffuse light that is directly reflected to the sensor. In this case, the scaling factor ( − )⁄ considers the real portion of the sky dome viewed by the target, with representing the slope angle between the inclined surface and the horizontal plane. The third term of the second member accounts for all the contributions diffused by the background. An exhaustive description of the model that includes all the steps to derive it can be found in Santini and Palombo (2019).
When the topographic correction is not selected, ( ) = ( ) and = 0. In this condition, the terms including Ad, compensate for each other, and a simpler model can be considered: If the correction of the adjacency effect is not selected, = . When no topographic and adjacency effects are selected, the model regresses to the simplified and most well-known model: where = + . Different models are used when different combinations of correction options are selected.

MODTRAN Simulation
The physical quantities of the different models are simulated based on the MODTRAN 4 model. Although MODTRAN 4 works only with homogeneous surfaces, its output allows one to discriminate between the radiance reflected by the target and the radiance diffused by the background and also provides a separate output for reflectance contributions from only the direct incident light. Based on this output, ImaAtCor keeps all the terms of the second member of equation (1) separate.
MODTRAN allows one to manage many parameters to simulate the necessary radiative quantities, based on the sensor and platform characteristics, the illumination and acquisition geometry, the surface elevation, and the atmospheric conditions and characteristics. MODTRAN can also simulate the quantities in a wide range of frequencies with a resolution as low as one cm -1 and provides several options for atmospheric scattering simulation, including multiple scattering models [16]. Depending on the correction selected, ImaAtCor operates a different set of MODTRAN simulations to estimate all the physical quantities of the Equations (1)-(3) ( , , , , , and other quantities needed for the estimation of the environmental contribution). The ImaAtCor GUI allows the user to manage all the MODTRAN input parameters needed to properly model the atmosphere and gives the user a choice between different operation modalities, as described in the section dedicated to the GUI.

Environmental Function
When the adjacency effect correction is selected, the is first estimated based on the simplified model of equation (3). Then, is computed on the basis of an environmental function ( ), which describes the efficiency of the atmospheric transmission of the signal reflected by a point located at a distance from the target. In polar coordinates, for a nearly nadir view, this relation between and ( ) can be expressed by the following equation from [17]: The environment function accounts for the probability that a photon within a distance from the target will be scattered to the sensor. Its derivative, ( ) , is known as the atmospheric point spread function. The ( ) depends on many parameters, including the phase functions of the molecules and aerosols, their vertical distribution, and their optical thickness [7,16,[17][18][19][20]. The derivation of the function ( ) is completely described in [4]. After obtaining , Equation (1) is used to derive the final reflectance . ImaAtCor offers different options for estimating ( ) and , with different precisions and processing times. The most accurate modality allows the construction and rigorous application of the environmental function described in [17]. This modality is highly time consuming and should be used only after all other atmospheric parameters have been defined to avoid unnecessary time expenditures. Once computed, can be stored for further applications. A much less time-consuming modality, but with a performance very similar to that of rigorous applications, can be obtained by applying a series of low pass filters to simulate the effect of the environmental function at various distances from the target.

Digital Elevation Model
When topographic correction is selected, the and angles of Equation (1) and the shaded areas are automatically computed by the tool on the basis of the DEM (digital elevation model) that is loaded by the user through the graphic interface. The computed parameters can then be saved in a single file (named associated DEM) and reloaded in subsequent applications to reduce the processing time. The DEM is also used as the input in MODTRAN simulations to account for height variations in the atmospheric column.

ImaAtCor GUI
ImaAtCor is provided with a user-friendly and exhaustive GUI, which allows the complete management of all data and implemented features. This GUI is organized with a main panel ( Figure  1) divided into sections for different functionalities. Each section allows the user to insert or select data and operate his or her choices through a series of dedicated menus and panels. This section is dedicated to a description of the GUI and its ability to manage the main features implemented in ImaAtCor. The following list briefly describes all the sections of the main panel: • "File": manages the input and output of files.
• "Conversion factors": allows the setting of conversion factors.
• "MODTRAN Models": manages the setting of atmospheric models for MODTRAN simulations. • "Sensor parameters": manages the setting of the instrumental parameters.
• "DEM Processing": manages the processing of the DEM for topographic correction.
• "Rrs processing": allows the computation of remote sensing reflectance for water bodies.
• "Template manager": manages templates (also for multiple runs in batch mode).
• "Correction levels": allows one to set the corrections to be performed and start the correction process. • "Parameter control table": shows the set parameters.

File
The "file" section allows the user to open the radiance file to be processed and select the reflectance file to save the result of the atmospheric correction. This tool provides a utility to automate the choice of the output file name through a series of suffixes corresponding to the different parameters considered during correction. This is particularly useful for multiple processing (see the dedicated section). Moreover, the automation of the naming process drastically reduces the risk of unintentionally overwriting the files and allows the user to always immediately obtain correct information about the content of a file starting from its name. The user is allowed to select the various suffixes to be used in the automated naming process through a choice menu ( Figure 2). The first suffix added to the base name is "_refl" (reflectance) and is mandatory. The other selectable suffixes relate to the correction levels applied to the image and to the values of the atmospheric and image acquisition parameters used for the correction: i) "_adj" refers to the correction of the adjacency effect; ii) "_topo" refers to the topographic correction; iii) "_rrs" refers to the remote sensing reflectance (used for the water's surface); iv) "_vis_xxx" refers to the values of the visibility; v) "_H2O_xxx" refers to the atmospheric column water content, etc. The choice menu with the complete set of suffixes to be selected appears when the F12 key is pressed.

Acquisition Parameters
The "acquisition parameters" section allows the setting of the geographic coordinates, the date and time, and the other acquisition parameters. The settable parameters are: • Date (dd/mm/yyyy).
The acquisition date can also be set using a date picker that appears when the calendar" button is pressed. If geographic information is present in the header associated with the image, the latitude and longitude values of the image center are automatically loaded and shown in their proper fields. Starting from the image latitude and longitude, the height is automatically computed using world DEM data stored in the "world_dem_geo_jpg.jp2" file.

Conversion Factors
This section allows one to set the scale factor to convert the radiance image into a MODTRAN radiance unit (µW/(cm 2 *nm*sr)). This form also allows one to set the conversion factor for output reflectance images. If the reflectance factor is set to 1, the output image will be a floating point; otherwise, it will be an unsigned integer with the reflectance value multiplied by the set reflectance factor. The "calculate factor" button calculates the radiance conversion factor starting from the input image radiance units, which must be set in the left column of the drop list menu (Figure 3).

MODTRAN Models
The "MODTRAN models" section allows the user to select the aerosol model and one of the six geographical-seasonal atmospheric MODTRAN models using the "atmosphere drop list menu". The VIS form allows one to set the visibility.

Sensor Parameters
The "sensor parameters" section allows the user to set the spectral responses and to define the look angles associated with each across track pixel. This section also includes a drop-down menu for the selection of a known sensor for specific configurations or data pre-elaboration.
When clicked, the "look angles" button opens a panel ( Figure 4) to configure the look angles in two modalities: i) Opening a file containing the look angles for each pixel of the scan line, and ii) setting the FOV (field of view) and the number of across track pixels to allow the system to compute the look angles. The computed look angle configuration is stored in the file selected by the user. If the "do not use look angles" box is flagged, the system associates the view angle corresponding to the inclination of the sensor axis during the acquisition ("vew zenith" value set on the principal panel) with each image pixel. The "responses" button opens a panel to define the spectral response functions ( Figure 5). The user can choose to open a file containing the spectral response functions in the ENVI (environment for visualizing images) spectral library format ("load from file") or to open a file containing two columns with the wavelength of the band centers and related FWHM (full width at half maximum). In the latter case, the spectral response function of each band is built as a Gaussian function centered at the band center and with a standard deviation of σ = √ . After selection, the spectral response functions of all bands are stored in the ENVI spectral format and plotted in the appropriate display ( Figure 5). The drop-down menu ( Figure 6) allows the selection of a specific sensor within a list. When a sensor is selected, the system automatically loads or computes the appropriate spectral responses and asks the user to select a file (generally provided with the image data) with the acquisition parameters to be loaded (i.e., MTL.txt for Landsat 7/8 multispectral or panchromatic images and .xml for Sentinel 2 images). The loaded parameter values are shown in the appropriate fields of the main interface panel (Figure 1). Within the selectable sensors, there are two modes referring to PRISMA: "PRISMA_ENVI", and "PRISMA_HDF5". When "PRISMA_ENVI" is selected, the user is asked to select a PRISMA image in the ENVI format and is allowed to select the original .he5 file to automatically load the acquisition parameters. The "PRISMA_HDF5" option allows one to directly use the PRISMA proprietary format .he5 through a dedicated interface, the "PRIMA reader tool" (Figure 7).
The PRISMA reader tool allows the user to convert the images from the .hdf5 to ENVI format. There are four options to convert the data: i) "Merge via SNR", which merges the VNIR and SWIR images by selecting the overlapping bands on the basis of the highest nominal SNR value; ii) "merge and sort", which merges the VNIR and SWIR images without operating any selection but sorting all the bands by wavelength order; iii) "merge only", which creates a merged image with the VNIR bands preceding the SWIR bands; and iv) unmerge, which saves the VNIR and SWIR images in two separate files in the ENVI format. The ENVI header file is built to include a series of ancillary information, such as wavelength and FWHM. By pressing the "save all ancillary data" button of the PRISMA reader tool, all ancillary data contained in the .he5 file are stored in separate files in ENVI format with the same name used in the original .he5 file. Finally, the tool also allows the user to operate spectral shift for the VNIR and SWIR separately to correct for a possible shift of the original data.

Settings up the MODTRAN Cards
In this section of the main graphic display, it is possible to manage all the MODTRAN cards to be used for the atmospheric simulations. Each button of the interface opens an interface for a specific card. In Figure 8, an example of a MODTRAN card setting interface is shown. By pressing the "cards reset to initial values" button ( Figure 1), all MODTRAN cards are reset to the values stored in the "card.ini" file, which is updated each time the main graphic display is closed.

DEM Processing
In the "DEM processing" section, the "associate DEM" button opens the "associate image-DEM" interface ( Figure 9). This interface allows one to select a DEM (georeferenced DEM), which is used to compute the other parameters needed for operating the topographic correction. The result is a multilayer image, referred to as an 'associated DEM', made up of five layers. The first layer constitutes the actual elevation model, while the other layers consist of the derived parameters: ii) sun angle; iii) shadow; iv) slope; and v) image mask. The 'associated DEM' is stored in a file named radiance file name + _dem and can be reloaded for further processing through the "choose associated DEM" button. The "associate image-DEM" interface offers three options to manage the shadow correction: i) "Calculate the shadow", which uses ImaAtCor to calculate the shaded areas based on the uploaded DEM; ii) "open the shadow file", which allows the selection of a file containing a mask of the shaded areas; and iii) "skip the shading correction", which skips the shadow correction.

Rrs Processing
The "Rrs processing" section allows one to calculate the remote sensing reflectance, which is the quantity usually used for water surfaces. This quantity excludes the contribution of the light reflected by the water surface and considers only the signal from the water column [4]. The software, in beta version, currently implements only an empirical estimate of the Rrs. This estimate is carried out by subtracting, from a reflectance image, the average value recorded in a certain range in the infrared region, where the Rrs is expected to be almost zero.
The "Rrs estimate" button gives access to the interface to estimate the Rrs value ( Figure 10). This interface allows one to load and display ("load RGB") a reflectance image. Using the left mouse button on the displayed image, it is possible to draw a ROI (region of interest) on a water area. At this point, the spectra of the ROI are plotted on the display, and the user is allowed to select the SWIR interval for the Rrs estimate using the mouse buttons. The converted Rrs image is then stored in the "Rrs output file" by pressing the "save" button.

Template Manager
The "template manager" allows one to save a complete set of run parameters in a template file; the "save template" and "restore template" buttons can be used to save and restore the template information. It is also possible to open several templates simultaneously to be processed in batch mode, which is useful both in the case of redundant applications on the same image and when a temporal series of images have to be processed.
ImaAtCor is also equipped with a template generator tool (Figure 11), which generates a series of templates that vary some parameters in a given interval. The tool opens when the F11 key is pressed. Template generation is initialized by a base template that has to be opened by the user. Then, the user is allowed to define the intervals and the steps of the parameters he or she wants to vary. The tool then generates a template for each combination of the set parameters, allowing the user to concurrently test a wide spectrum of correction parameters in batch mode.

Correction Levels
This section allows one to select the correction to be performed and starts the correction process. Any combination of adjacency and topographic correction can be selected, ranging from a base correction (if none of them are selected) to complete correction that includes both adjacency and topographic corrections and corresponds to the model in Equation (1).
When the adjacency correction is selected, an appropriate interface ( Figure 12) allows the user to choose between two different modes for applying the environmental function (see 2.3), which differ in their rigor and processing time: i) "Rigorously apply Equation (1)", which is very expensive in terms of calculation time, and ii) "perform a series of low-pass filters instead of computing the integral of the equation", which significantly reduces the processing time with a small loss of precision. In both methods, the environment function is computed following [17], as explained in [4]. The interface also allows one to set the spatial resolution, which is automatically set when a geocoded radiance image is open. This value is used to compute the environmental function.

Control Table
In the parameter control table (Figure 13), the user can check the input of the correction parameters. Red refers to mandatory parameters, and yellow refers to optional parameters. Red and yellow are replaced by white after the parameters have been set.

ImaAtCor Application Results
This section illustrates some noteworthy results of the ImaAtCor application relating to topographic, adjacency, and shadow corrections. Figure 14 presents a subset of a normalized spectral distance image [4] between two LANDSAT-8 images acquired from the arid region of Parinacota (Andes, Chile) in the same location but during different seasons (different illumination geometries) and corrected excluding (left) and including (right) topographic correction. The image shows remarkable flattening when topographic correction is included. To demonstrate adjacency correction, Figure 15 depicts a subset of a LANDSAT-8 image of Pertusillo Lake (Lucanian Apennines, Italy), a small inland water body surrounded by flourishing vegetation. In the same figure, the graph on the left shows the mean reflectance spectrum of the ROI on the lake (yellow rectangle), obtained while excluding (dashed line) and including (solid line) adjacency correction in the correction procedure. The graph clearly shows the presence of a sharp adjacency effect due to the near-infrared peak of vegetation on band 5, in which the water shows almost no signal due to water column absorption [5]. After the adjacency correction, the impact of the surrounding vegetation is greatly reduced. As the Rrs estimate was not included in the correction, a residual value in the infrared range is still present due to the water surface reflectance. The Pertusillo image (centered at 40°10'0"N;16°0'0"E and 40°17'0"N;15°57'0"E) was acquired on 2019/08/24 at 9:35 UTC.
In Figure 16, an airborne CASI reflectance image acquired from an industrial area in the province of Reggio Calabria (Calabria, Italy) is shown before (Figure16a) and after (Figure16b) the application of topographic and shadow corrections. The results are remarkable, as the correction allowed almost complete recovery of the signal in the shade. This is confirmed by the graph in Figure16c, which shows the spectra corresponding to an exposed and shaded portion of a paved road before and after shadow correction.
(a) (b) (c) Figure 16. CASI (compact airborne spectrographic imager) reflectance images before (a) and after (b) the application of shadow correction. The ROI (region of interest) of the zoomed area relates to the shadow (red) and exposed (blue) areas. The graph shows the mean spectra of the ROIs before and after the application of topographic corrections (c).
For a more detailed description of the above results, refer to [4]. As a further demonstration of ImaAtcor's potential, Figure 17 presents the results of applying adjacency correction on a soil target. Figure 17a depicts an agricultural area in the province of Metaponto (Basilicata, Italy). The graph of Figure 17b shows the reflectance signal at 871 nm, corresponding to very high vegetation reflectance values, along the transept (red arrow) drawn on Figure 17a before and after correction for adjacency effects. A 10% increase in the contrast index was obtained (calculated as (max-min) / (max + min)); as expected, the intensity increased in the vegetated areas but decreased in the bare soil. The correction was applied to a non-georeferenced image to avoid any kind of distortion that would have made the image less sharp, partially hiding the effects of the correction. For this reason, the abscissa of the graph shows the number of pixels from the beginning of the transept. The Metaponto image (centered at 40°23'00"N;16°50'28"E) was acquired on 2013/07/16 at 14:59:30 UTC at 3200 m above sea level, corresponding to a spatial resolution of about 1.50 m.
To determine the processing times, we carried out a series of typical runs on a standard 260band PRISMA image on a PC equipped with an Intel(R) Core(TM) i7-7700HQ CPU @2.80 GHz. The results are shown in Table 1. The third row of the table refers to a case in which the correction of adjacency effects is included but a previously calculated image is reused.

Conclusions
This paper presents the ImaAtCor tool, including its various models and GUI. ImaAtCor is a versatile tool for the correction of remotely sensed data for atmospheric and topographic effects. The implemented models are physically based and take into account second-order corrections, such as adjacency effects. The software's exhaustive and friendly GUI allows the user to manage the correction process, including MODTRAN simulations, to estimate the necessary atmospheric quantities. This tool includes several utilities for the pre-processing of PRISMA data and provides an instrument for the management of multiple data in batch mode, which is very practical when temporal series have to be pre-processed. The implemented methods demonstrate great potential, as shown in [4]. ImaAtCor is used on a daily basis to correct and enhance the quality of images in research contexts and experiences continuous updates. In the near future, the implementation of a physically based method for the derivation of the Rrs is planned.

Software Distribution
This software is currently distributed on demand. It is possible to request the beta version of ImaAtCor for evaluation purposes directly from the authors at angelo.palombo@cnr.it or federico.santini@cnr.it. The tool was implemented in a Windows environment and requires an IDL virtual machine 8 or higher and the installation of MODTRAN 4.
Author Contributions: Investigation: F.S. and A.P., methodology: F.S. and A.P., software: F.S. and A.P., All Authors contributed to the preparation, review and editing of this manuscript. All authors have read and agreed to the published version of the manuscript.