17.2. Technical Description

17.2.1.1. Traditional Photon Mapping

The new approach that was identified is a modeling technique called "photon mapping" [Jenson 2001]. The photon mapping approach is a hybrid of traditional forward and backward Monte-Carlo ray tracing techniques. In this two-pass method, source photons are shot from a source into the scene using forward ray tracing during the first pass and then collected using a backward ray tracing during the second pass (Figure 17-1). In the original presentation, photon mapping is described as a monochromatic approach, however various strategies have been introduced to handle spectral sources using this methodology. This technique has been demonstrated to be applicable to traditional solid geometry reflective illumination and scattering and absorption by participating mediums.

During the first pass, a photon performs a pseudo-random walk through the scene based on the optical properties of the surfaces or mediums that it encounters. At the location of each interaction, information regarding the interaction event is stored into a fast, 3D data structure (e.g. a kd-tree) that can be queried during the second pass. The modeled photon is followed until it is absorbed by some element within the scene, and the process is repeated for some number of modeled photons that originate at the source. During photon forward tracing pass, the characteristics of the source can be employed including the directional characteristics and the spatial, spectral and temporal distribution of the source photons.

During the second pass, rays are shot from the image plane into the scene until a scene element is intersected. At that intersection point an estimate of the incident flux must be computed. This is accomplished by querying the 3D data structure that was filled in the first pass to find a set of incident photons arriving at the surface (or volume element in the case of a participating medium). Each of these photons is then then redirected back towards the sensor using the applicable surface reflectance model or volume scattering model using the photon's specific incident and excitant geometry (Figure 17-1).

From an absolute radiometry perspective, each modeled photon is equivalent to some number of absolute photons emitted from the source. Depending on the angular extent of the beam, the field-of-view of the sensor and the spatial detail of the scene approximately 250,000 to 1,000,000 photons can be utilized in the photon mapping process to create reliable statistics for topographic returns. For volume scattering returns, this number may need to be higher depending on the absorption and scattering coefficients and the spatial extent of the participating volumes.

17.2.1.2. Photon Mapping Enhancements

In order to adapt traditional photon mapping for the purpose of active laser radar applications, a set of enhancements was made to the basic treatment. The primary addition was the addition of a field to track each photon's total travel time, which allows the travel time accrued during all of the multiple bounce/scattering events to be recorded and used for time gating. Future modifications may include tracking the polarization state, the phase distributions and the spectral structure of the photons represented by each modeled photon stored in the photon map. At this time, all interactions within the model are elastic in nature, but the incorporation of inelastic collisions can be handled within this treatment by imposing a combination of wavelength and temporal shifts during each interaction.

17.2.2.1. Forward Propagation

The photon shooting function leverages the generic ray tracing support that already existed within the DIRSIG model. The ray tracer interacts with scene elements that have material specific properties. Each material has a set of surface optical properties and an optional set of bulk or medium properties. The surface properties include a spectral reflectance and/or emissivity property. The currently supported reflectance (BRDF) models include importance based sampling functions to support the forward and backward Monte-Carlo ray tracing. The bulk properties include spectral extinction, absorption and scattering coefficient models. When volume scattering is modeled, a phase function object can also be configured to describe the directional nature of the scattering. The phase function objects also support importance based sampling functions to support the forward and backward Monte-Carlo ray tracing. Although these optical descriptions existed prior to the addition of the LIDAR capabilities, their interfaces were enhanced during this effort to facilitate an efficient implementation of the photon mapping subsystem.

17.2.2.2. Backward Propagation

The next major component of the model to be implemented was the collection of photons from the photon map and forward propagation into the imaging system. The modeling of the radiometric returns from a scene element is handled by a class of objects in DIRSIG referred to as radiometry solvers. A radiometry solver encapsulates an approach to predicting the energy reflected, scattered and emitted by a surface or volume. For example, there already existed a specific solver for predicting the passive returns "from hard" surfaces and another for a volume of gas (e.g. plume). These solvers utilize the material specific surface and bulk optical properties to predict their results. One or more radiometry solvers can be assigned to each element in the scene. To support the active LIDAR returns, a new radiometry solver was created to compute the returns from a scene surface or volume by using the optical properties and the photon map to estimate the number of incident photons at the elements point in space. Unlike the existing passive radiometry solvers that would place the final result in a time-independent result object, the LIDAR-specific radiometry solver places the result into a time-gated result object. The time gating is based on a user-defined signal gate consisting of a start, stop and delta time.

The final component was the implementation of a "capture method" object within the system that directs how outputs for the detectors on the focal plane are computed. The LIDAR specific capture method developed under this effort backward propagates rays from each detector element into the scene where it intersects scene elements, runs the appropriate radiometry solvers (passive and active), forward propagates the energy to the focal plane and then writes the arriving photons counts to the output file.

17.2.2.3. Instrument and Platform Modeling

Most operational laser radar instruments are flown on aircraft and utilize some method of aircraft relative scanning to increase the spatial coverage of the system. The changes in viewing geometry during the scanning process and the location, orientation and stability of the instrument platform can affect the final data products. For example, the ability to resolve a specific object in a topographic data product derived from a dataset might be dependant on the angle from which it is illuminated. The overall accuracy of a derived topographic product might depend on the overall stability of the platform and knowledge of the platforms position. The DIRSIG model has a flexible platform model that allows the platform to be positioned and oriented as a function of time. Furthermore, the instrument can be pointed with respect to the platform either statically or dynamically using one of the available instrument mount objects. These mount objects can support temporal scanning including basic sinusoidal across-track scanning as a function of a user-defined scan rate.

17.2.2.4. Generalized Modeling Process

With the basic components of the model now described, the overall modeling process can now be summarized. A modeling run consists of the user specifying the scene to be modeled, the instrument and instrument mount description, the source description, the platform positioning data, and a set of tasks that describe time windows over which data is to be generated. The data generation process begins by walking through each user-define task according to a step time that is usually driven by the source pulse rate. During each time step, the platform and instrument mount positions and orientations are computed for the current time, the source is fired which fills the photon map, the focal plane is captured which collects the photons and propagates them to the sensor, and the focal plane reaching photon counts are written to the output file. This cycle repeats for each time in the task window and for each task in the list. The final product of the tool is a 3D cube consisting of photon counts as a function of two horizontal spatial dimensions and one temporal dimension. This data cube is ingested by external focal plane models to handle instrument specific detection schemes (e.g. linear versus Gieger mode), noise sources, etc. Further external processing of the resulting modeled instrument outputs can be used to create final data products (e.g. topographic maps, gas detection maps, etc.).

17.2.3.1. Atmospheric Optical Structure

At this time the atmosphere is assumed to be spatial uniform as a function of position and altitude, which is not accurate for many real applications. The horizontal and vertical structure of the atmosphere results in different absorption and scattering characteristics as a function of location. Ideally, robust atmospheric optical models like MODTRAN and FASCODE would drive both the extinction and scattering optical properties of the atmosphere. Both of these models are currently horizontally homogeneous, however, external methods can be used to horizontally combine the vertical structures of these models to create a pseudo-3D atmospheric model. At this time, the extinction coefficients used by the DIRSIG model are extracted from the existing MODTRAN and FASCODE derived tables. However, extraction of the vertically structured scattering coefficients and phase functions from MODTRAN (FASCODE does not support scattering) would require custom modifications to the MODTRAN code. A request for these modifications to the MODTRAN distribution has been requested. For the time being, the user can create their own scene elements with user-defined extinction and scattering properties to replace the atmospheric if this level of control is critical.

17.2.3.2. Atmospheric Returns

Under most conditions, the extinction and backscatter coefficients of the atmosphere are extremely small, which means that the probability of absorption and scattering events within the atmosphere is very low. For example, the scattering coefficient for a dry atmosphere might be 10e-5 [1/m], which means that you would need to shoot 10e+5 photons into a 1-meter long box of atmosphere to observe a single scattering event. Many systems are attempting to resolve vertical resolutions of a fraction of a meter and from an altitude of several thousand meters, which implies that you would need to model 10e+10 photons within each spatial detector element in order to get a single scattering event within each numerical contribution element. To achieve robust statistics, this number would be ideally several orders larger. Therefore, the numerical model approach using photon mapping is not considered a tractable for this part of the problem.

To use the numerical approach utilized by the photon mapping technique, the number of photons that would need to be shot into the atmosphere to create a statistically accurate representation of the scattering events would be many of orders of magnitude larger than the number of photons needed to model the topographic returns. To avoid the problems of predicting the atmospheric returns numerically, the atmospheric returns from the model are currently modeled analytically using the formulation proposed by Measures ([Measures 1984]):

Equation 17-1.

where P is the total scattered laser power received in watts, is the peak wavelength of the laser in meters, R is the range of the contribution element in meters, is the average power of the laser pulse in watts, c is the speed of light in meters per second, is the pulse length in seconds, is the receiver optics transmission at the laser wavelength, is the geometrical form factor at range R, is the atmospheric backscatter coefficient at the laser wavelength and range R in inverse meters and is the atmospheric total extinction coefficient in inverse meters. For most of the systems modeled to date, the atmospherically scattered photon returns amount to only a few photons accumulated over the entire path length, which is far below the detection level of the modeled systems. However, if the tools is used to model significantly longer path length or an optically thicker atmosphere (containing fog), then these numbers will grow to be large enough for consideration by the detection model. The current output of the model is integer photon counts for each time resolution bin, which results in truncation of these photon counts. For these situations, the DIRISG model should be modified to output fractional photon arrival counts so that the detector model can properly accumulate them.

17.2.3.3. Temporal Pulse Shape

The temporal shape of the source pulse is model parametrically at this point in time. This saves a significant amount of run-time and memory resources by allowing the model to track pseudo-photons that represent temporal bundles of actual photons. If this approach were not employed, the number of model photons would have to be dramatically increased to correctly sample the temporal envelope of the pulse. The parametric temporal shape approach will break down when inelastic collisions are introduced which might skew the temporal shape of the pulse reflected or scattered by such an material. In addition, the user-defined signal gate should be setup to provide accurate sampling of the temporal pulse shape since it is assumed than many detector systems will detect or trigger on the rising slope of the pulse. Course temporal sampling of the pulse can disable the uncertainty in the probability of detection within an external focal plane model unintentionally because the photon counts can change from far below the detection threshold in one time bin to far above the detection threshold in the next time bin. Ultimately, the temporal resolution of the signal gate defined by the user can limit the height accuracy of a derived topographic map.