NOTE: Typically the Input image has no map coordinates
associated with it, since this is the image to be reprojected.
Because of a glitch in Imagine, it is necessary to go through
the following steps for images that have no map coordinates
associated with them before the start point is in the viewer.
Open the input image in a viewer (section IV, steps 1–2
Open the ImageInfo for that image (section IV, step 3
Select “Edit/Change Map Model.” The Change Map
Model window will appear. Simply click on OK button
once this has been done.
Accept changes when prompted.
Close the viewer.
Why do this? For some reason, Imagine has a
disconnect between the map info that is used in the
viewer when an image is opened that has no geometric
information and the map information that is passed to
the ITPFind program and to the GCP editor of Imagine
(further down in this SOP). By opening the Edit/
Change Map model dialog, this forces Imagine to make
geometric properties consistent across modules.
Pixel Center to Pixel Center Distance:(x,y):
The distance in a COMMON REFERENCE SYSTEM
between the center of adjacent pixel centers in this image.
When ITPFind runs, this value will be compared to the
parallel measure in the complementary reference or input
image to determine how far the program should jump between
ITPs. It must be in the same units for both the reference
and the input images in a given run of the program, typically
in units of the reference image. Since the geometric-reference
image is a LANDSAT TM image with pixel size 25 by 25
m, its pixel center to pixel center distance should be listed as
25, 25 m. However, assuming that the input image has been
ordered from EDC with 28.5 m pixels, its pixel center to pixel
center distance should be 28.5 m, even if no map information
is associated with the image. Note that if the input image
coordinate system is unknown (a likely case, given that ITPs
are being sought), the “pixel size” of the input image may be
reported in the ImageInfo of ERDAS Imagine as 1.0, 1.0 – but
the Pixel center to pixel center distance:(x,y) is still 28.5, 28.5.
The rotation, in degrees clockwise, of this image relative
to a stable reference direction. It typically is easiest to set
the reference image rotation to zero, and set the input image
relative to the reference image. To determine whether the
rotation of the input image relative to the reference image is
positive or negative, use the following approach:
1. Imagine an arrow pointing north on the ground in the
2. Imagine an arrow pointing north on the ground in the
3. If the head of the input image arrow is to the left of the
reference arrow, the input image rotation is negative. If
the head of the input arrow is to the right of the reference
arrow, the input image rotation is positive.
Layer to Use:
The layer number of the image to use for the analysis.
Only one layer is used. In all likelihood, this should be the
same layer number as will be used in the complementary
reference or input image. Layers begin with 1. It is best to
use a layer of the image with high contrast, such as the near
infrared bands. For TM images, band 5 is often good.
If the image file contains unwanted background, give the
value of the background here. This is often zero. The programs
will recognize background areas and will skip past them,
improving performance. Moreover, this method can be used to
“work around” popcorn-type clouds or other drastic changes
between images. A quick unsupervised classification of the
cloudy image can be used to identify the bulk of the clouded
areas and the shadows, and all of these can be set to the Ignore
value. The program will use only areas between the clouds
for calculation of correlation. Note that if this route is chosen,
the window size parameter (see below) likely SHOULD be
increased to ensure that an adequate number of pixels remain
for calculation of correlation when the cloud pixels are
Mask Below, Mask Above (Both Optional)
As of version 2.1, a more useful approach to screening
out unwanted areas was added by inclusion of the Mask Below
and Mask Above keywords for an image. These are optional
– if not used, do not include the keyword on the line (i.e.
simply omit the entire line, rather than retaining the keyword
“Mask below” and setting it blank). The number attached to
either mask below or mask above is the value below or above
which pixels will not be used for correlation matching. This
was designed for use in cloudy areas, where a mask above
value eliminates much of the cloud area and the mask below
value eliminates much of the deep cloud shadow. The values
must be chosen by the user, so it is best to just look at your
image in your favorite image processing system and get a
rough guess for the digital number (DN) values of the clouds
and shadows. As with the ignore keyword, it is best to increase
the window size you are using if you think that a fair portion
of each window will be masked out with this feature.
SOP 2. Preprocessing Landsat Imagery
1. Parameter files are pointed to in the “File of available
parameter files” (see above). Each run of the software
requires a single-parameter file. An example of the
required information is given below, followed by detailed
descriptions of each line.
Original Parameters File for use in ITPFind program
Window Size: 150,150
Window Spacing: 400, 400
Number of Iterations: 1
Pixel Aggregations: 1
Search Neighborhoods: 5
Threshold min. steepness: .35
Zoom Factors: 2
Maximum move: 0
The size, in pixel counts of the reference image, of the
window extracted for matching. It is in x,y format. Note that
this dimension is not in the units of any coordinate system, but
referenced by the number of reference pixels. For example,
if the reference is a TM image with pixels 28.5 m on a side,
then setting the window size to 100,100 means that a window
of 100 by 100 pixels will be extracted, equivalent to an area
2,850 by 2,850 m. The equivalent sized window in the input
image is calculated by the software based on the
to pixel center distance:(x,y)
given for the reference image
and the input images. If they are not the same, the program
will resample the input image to match the
Pixel center to
pixel center distance:(x,y)
of the reference image.
The approximate desired grid spacing between ITPs,
in pixel counts of the reference image, between ITPs.
Again, the relation between the
Pixel center to pixel
of the reference and input images will
determine the jump in pixels of the input image. Note that
the rotation of the input image also will be calculated in the
jumping between points. Thus, if the user is unsure of the
rotation value of the input image, it is safer to give a smaller
window-spacing value, since the error in the rotation value
will be multiplied over a shorter distance.
Number of Iterations:
The number of times that the program will refine the ITPs
coordinate pair at each point. The next five parameters must
each have this many entries, each separated by a comma. For
Number of Iterations
is set to 2, then each line
of the next five must have the keyword parameter followed by
“: <n1>, <n2>”, where <n1> and <n2> are the values for that
Typically, this can be set to 1. If set to 2 or higher, the
program engages in an iterative process. Using the first
entry in each of the next 5 parameters, the program finds an
approximate match for the ITP pair. Using that approximate
match derived in the first iteration, the program uses the
second entry (the number after the first comma) to find a better
match for the ITP pair. Typically, this only makes sense if the
value is lower for the second iteration
than for the first iteration.
The number of reference pixels by which to aggregate
both images before floating the input image over the reference
image. Note that this aggregation occurs AFTER the input
image is resampled to match the grain size of the reference
If set to 2.0, for example, then the reference image
is first aggregated such that a 2 by 2 set of original pixels
becomes 1 pixel in an aggregated image. That image—with
its doubled pixel center to pixel center distance—is used to
develop correlation surfaces and find an ITP. Because of the
aggregation, this ITP will have a lower precision than if the
ITP were calculated on the images at native resolution, but this
aggregation allows a larger area, with more potential for strong
spatial pattern, to be used to calculate correlations. Used in
conjunction with two or more iterations (set with the Number
of iterations parameter above), this can efficiently hone in on
ITPs in difficult-to-match image pairs.
The number of pixels to float the input image is relative
to the reference image. This value is in units of the aggregated
pixels—i.e. if the
for this iteration were
2, then a search neighborhood of 10 for this iteration actually
would mean that the input image has floated an equivalent
of 20 original pixels in all directions relative to the reference
This parameter greatly influences performance of the
software, since increases in the search neighborhood result in
squared increases in the number of comparisons necessary. It
is thus desirable to keep this number as low as possible. Values
between 5 and 7 are optimal, although values as high as 10 or
12 may be necessary.
6 Protocol for Landsat-Based Monitoring of Landscaped Dynamics at North Coast and Cascades Network Parks
Threshold Minimum Steepness:
This value is the threshold against which a potential
peak in correlation is tested, based on the relation between the
peak of correlation and the plane at the “base” of the peak. If
the observed value is below this
Threshold min. steepness
value, the point will be thrown out. A value of 0.35 has been
robust across most image situations. If it seems that erroneous
ITPs are being accepted too often, increasing this parameter
value may help—however, generally it is better to attempt
altering other parameters first, i.e. the window size or the
The portion of the originally extracted window
(determined by the window size) is used for the determination
of covariance. A value of 1 means that the entire extracted
window is used to determine covariance. A value of 2 means
that the window used for covariance calculations will be
one-half the size of the original image in both the X and the
Y dimensions; i.e. it will be one-quarter the area. The smaller
window is extracted from the center of the originally extracted
Consider the following example: The window size is set
to [400,400], number of iterations set to 2, pixel aggregations
set to [2,1], and zoom factors set to [1,2]. In iteration 1, the
400 by 400 window is read from the file, and then aggregated
in 2 by 2 blocks to create a 200 by 200 image that captures
the spatial patterns across the entire area represented in the
original [400,400] image. The spatial matching is run on
that image and a preliminary ITP pair located. In the second
iteration, the zoom factor of 2 takes the original 400 by 400
image and extracts the 200 by 200 area around the preliminary
ITP pair found in iteration 1, and this image (unaggregated,
because the pixel aggregation is set to 1 for this iteration) is
used to derive a more refined ITP pair.
If this value is non-zero, every time the program cannot
find a valid covariance peak at a given point, it will try again
by searching nearby rather than skipping the point altogether.
“Nearby” is defined by the next parameter, the
The pattern of searching nearby is circular; i.e. it will attempt
to find a point just to the “right” of the first attempt, then
slightly right and down if no match was found there, then
down if no match is found, etc., until a complete circle has
been traversed. If the
is 2, this process will
continue on a second circle with greater radius.
In most cases, it is more efficient to set this value to zero,
since factors that cause difficulty locating peaks at a given
location (for example, a bank of clouds in one image) are more
likely to exist nearby than far away. However, under certain
circumstances the user may require that ITPs be located as
close to full-grid pattern as possible, and this option allows
This value determines how far away from the original
point the program should search for covariance values if it
cannot find a covariance peak on an initial try. The value is
a proportion of the window-size value. If the window size is
[400, 400] and the
is 0.5, then the program will
extract new windows at a radius of 200 reference pixels away
from the original point tried.
Some tie points will still be odd, despite the various
filters that occur during the search for each individual point.
By setting this parameter, the program calculates a simple
first-order transformation after all points have been found,
and iteratively removes any points whose removal results in
an improvement in the overall root-mean-square (RMS) error
of the solution. The value of this parameter is the threshold
for determining whether improvement has occurred, and
corresponds to the proportional improvement in RMS error.
Thus, setting it to 0.025 means that points will be removed
when their removal improves the RMS error by 2.5 percent or
more. This tends to remove the points that clearly are errors
(i.e. points over water, etc.). This feature was added with
At a minimum, the squeet_master.txt file must exist, and
it must point to a squeet_images.txt file and a squeet_params.
txt file. The squeet_images.txt file must contain information
for at least two images, one to be used as reference image
and one to be used as the input image. The squeet_params.
txt file must point to at least one parameter file with all of the
necessary parameter fields filled in.
Running the ITPFind Module
To find a good start point of the input and reference
images, load each image into an Imagine Viewer (follow steps
1 and 2 of section IV above). In the reference image, use the
side bars on the viewer to navigate to approximately the center
of the image, and visually locate a feature on the landscape
that is recognizable and likely to change little over the time
period of the two images. In the viewer, select the crosshair
icon (it looks like a plus sign, and is just to the left of the
hammer symbol) to start the coordinate-crosshair tool.
SOP 2. Preprocessing Landsat Imagery
The X and Y coordinates of the crosshair are shown in
the X and Y fields (fig. 14). Grab the center of the crosshairs
by clicking with the mouse point, and drag it to the feature
desired. In the viewer for the input image, use the sidebars to
navigate to approximately the same area of the image, locate
the landscape feature that was identified in the reference
image, start a coordinate-crosshair tool in the input image,
and click and drag its center point to approximately the same
location on the landscape as the crosshairs in the reference
Record the X and Y position of the respective crosshairs
in the squeet_images.txt file for the reference and the input
Once pointer files, image files, and parameter files have
been defined, drag the ITPFind icon over the IDL-VM icon to
start the program. It will query immediately for the location
of the squeet_master.txt file that was created above. Navigate
in the file-finding dialog box to the appropriate directory and
select the file.
An action window will pop up. The list of available
image files will be on the left, the list of available parameter
files on the right. On the bottom is an empty table. This
table is where reference and input images will be paired and
associated with a parameter file.
Click on an available image that is the reference-
geometric image. Its name will appear in the ‘Currently
Selected’ box. Then click in the first box (first row) of the
reference column (the left-hand column) in the bottom table.
For the reference image in ITPFind, use the 25-m
geometric-reference image that has been clipped to the
study area (section IV, step 12 above). This will ensure that
tie points are located only within the study area where terrain
information in the form of the DEM also is available.
Pick the input image from the same available image list
in the same manner and place it in the middle column below,
and finally pick a parameter file and place in the parameter
column below. At this point, the ITP program can be run by
clicking on the ‘Submit for processing’ button.
Figure 1. Coordinate-crosshair window.
Alternatively, additional sets of reference/input image and
parameter files can be placed in the second row of the bottom
table. After the ITP program has found ITPs for the first set, it
will move to the second set, third set, in order.
While the program is running, it will display two
windows: one shows the growing grid of ITPs, the second
shows the image chunks and correlation surfaces for the points
being processed currently.
The program will output five text files: one file for X
coordinates of the input image, one for Y coordinates of the
input image, one for X coordinates of the reference image, one
for Y coordinates of the reference image, and one readme file.
All files are placed in the directory of the INPUT IMAGE.
They will be named with an identical root based on the
image codenames used to find the ITPs and the time at
which the program commences. The files with the X and
Y coordinates will be imported into the ground-control-point
(GCP) editor of Imagine below (see step 4.b. below).
These files are only the coordinates of the ITPs. From this
point, Imagine must be used to conduct the terrain correction.
Conducting Terrain Correction
1. Load the geometric-reference image that has been clipped
to the study area.
a. The image was created in section IV, step 12 above.
b. To load the image into the viewer, see section IV,
steps 1 and 2 above.
2. In the main icon panel of Imagine, select the viewer
button (fig. 15).
a. Viewer #2 will come up.
i. Open the Input image in Viewer #2.
(1) Again, click on the folder icon of the viewer
and navigate to the appropriate folder. This
assumes that the input image has been
imported into Imagine format. If not, it
is necessary to import this image into the
Imagine format (section III above).
Figure 1. A portion of Imagine’s main icon window,
with the button that starts Viewers circled.
Protocol for Landsat-Based Monitoring of Landscaped Dynamics at North Coast and Cascades Network Parks
3. From within Viewer #2 (with the input image), select
a. An action window named “Set Geometric Model”
will pop up, with a list of options. Another action
window labeled “Geo Correction Tools” will pop up
(see fig. 8 above).
b. In the Set Geometric Model window, select Landsat
and click on the OK button.
i. The “Landsat Model Properties” window will
pop up (fig. 16).
c. For sensor type, ensure that TM is selected, and that
the Landsat number: field is correct.
i. The information on which Landsat sensor
acquired an image can be found in the header
file (a file ending with .h1) that came with the
original Landsat image from EDC (see SOP 1
d. In the “Elevation File: (*.img)” field, use the folder
icon to navigate to the DEM that has been clipped to
the study area (this step is discussed in section V, step
e. Ensure that the elevation units are correct.
i. This information should be found in the metadata
associated with the original DEM file.
ii. If there is no such information on the vertical
units of the DEM, it can be inferred in many
(1) Open the DEM in a viewer (see section IV
(2) Open the ImageInfo for that viewer (see
section IV step 3).
(3) In the Statistics info section of the
ImageInfo window, note the value of
the “Max:” field. This value should
approximately represent the maximum
elevation in your study area. Since units of
feet and meters differ by a factor of three,
this maximum value should indicate which
unit is used for the DEM. Other vertical
units (decimeters) have sometimes been
used; if neither feet nor meters seem likely,
then it will be necessary to track down
the source of the DEM and determine the
f. Leave the rest of the fields in the default position.
4. Setup the coordinates for the terrain correction by
importing the image tie-point text files that were created
by ITPFind above.
a. Click on the “Projection” tab of the Model properties
i. Click on the button labeled “Set Projection from
(1) The “GCP tool reference setup” window
will pop up, with the “Existing viewer” line
selected with a circular radio button. This
is the desired choice, so simply click on the
(a) The “Viewer selection instructions”
window will pop up, instructing user to
left-click the mouse in the viewer with
the reference image. Place the mouse
in Viewer #1 and left-click the mouse
(b) The Reference Map information
window will pop up, with a list of the
geometric properties of the reference
image. Click on the OK button.
(c) Imagine will then rearrange windows
and add windows.
b. All of the steps below occur in the “GCP Tool”
i. First, turn off the automatic calculation and
display mode by clicking on the “Toggle fully
automatic GCP editing mode” button on the left-
hand side of the window (fig. 17).
Figure 16. Landsat model properties page, used to reproject
images acquired by the Landsat sensors.
SOP 2. Preprocessing Landsat Imagery
(1) With the mouse over the title of the “X
Input” column, left-click on the mouse to
select that column. The first cell should be
highlighted and the shaded title box of the
column inverted, indicating that the column
has been selected.
(2) Right click from the same position over the
title box to see the column options window.
From the list of options, select “Import.”
(a) The “Import column data” window will
appear. Navigate to the directory where
the Input image used in the ITP find
program above is stored. This is where
ITPFind puts tie-point files.
(i) Select the file whose filename ends
with “…_xinput.txt” and click on
the OK button.
(1) Note: The first part of the
file name corresponds to the
image code names that were
included in the available image
text files used as input to the
ITP program, the second part
corresponds to the time stamp
when the ITP program was
(ii) The X coordinates of the input
image tie points should appear in
the X-input column.
(3) Repeat this process for the Y input column,
this time selecting the filename ending with
(4) Repeat for the X ref and Y ref columns as
(5) When all the tie points have been located,
the elevation values in the Z ref column
should appear. If not, click on the “Z” icon
in the GCP Tool window to load them. If no
values appear, return to the Landsat Model
properties dialog (fig. 16) and ensure that
the DEM has been identified properly, and
then click on the Apply button again.
(a) If no values appear after this step,
open new viewers and ensure that the
coordinates of the DEM correspond to
the same area as the image.
(i) Open a new viewer.
(ii) Load the reference image.
(iii) Load the DEM, but when loading
the image in the “Select Layer to
Add” dialog (fig. 18), click on the
“Raster Options” tab.
(iv) Uncheck the “Clear display” box.
Now the DEM will load into the
same viewer as the reference
(v) Once the DEM has loaded, right-
click in the viewer and select
“Fit image to window” from the
(vi) The viewer should redraw. The
image should appear first, then
the DEM should appear over
it. If the two are in entirely
different portions of the viewer,
this indicates that coordinates or
projections have been in error.
Double-check all resampling and
reprojecting steps above.
Figure 1. GCP tool, with the automatic GCP mode button
Figure 1. Window where files are chosen to add to a Viewer.
0 Protocol for Landsat-Based Monitoring of Landscaped Dynamics at North Coast and Cascades Network Parks
(6) Calculate the solution for the reprojection.
(a) Hit the sigma symbol on the left-hand
side of the GCP Tool window.
(i) The columns X residual, Y
residual, RMS error, and Contrib.
all will be populated. Additionally,
the field in the center top of the
GCP tool will display the “Control
point error” for X, Y, and total. This
is the total RMS error. See Richards
(1993) for more information on the
calculation of the RMS error. These
will be in units of the input image.
If the input image is without map
information and has pixels of size
1.0, for example, then a total error
of 0.5 would indicate a RMS error
of half a pixel.
(7) Screen out any tie points that are outside the
edge of the DEM.
(a) With the mouse over any row in the
GCP tool, right click to get the row-
selection tool. Select “Criteria” to bring
up the “Selection Criteria” window
(b) Using the mouse and the left mouse
button, select from the Columns: field
the “Z Ref.” row.
(i) Imagine will place $“Z Ref.” in the
“Criteria:” field at the bottom of the
window, and will continue to build
the criterion expression as more
values are selected below.
(c) Select from the “Compares” field the
(d) Select from the keypad the number 0.
Figure 1. Selection criteria window for the GCP tool.
(i) Ensure that the final expression is
$“Z Ref.” == 0.
(e) Click on the Select button. All of the
rows in the GCP tool where Z Ref. is
zero will be highlighted in yellow.
(f) Ideally, there should be none of these
rows. But in case there are, place the
cursor over any row number in the left-
hand side of the GCP tool and right-
click. Select “Delete selection.” This
will eliminate any points where the
elevation is zero.
(i) Why? Elevations of zero would
otherwise be incorporated in the
general solution relating the images
to each other, and false zeros will
cause large distortions in the
(8) Recalculate the solution for the reprojection.
(a) This must be done twice, once before
and once after taking out points with
zero elevation, or else Imagine will
(9) The control-point error should be well below
1 pixel in size, ideally less than one-half of
a pixel. If it is significantly greater than this,
examine the control points in the table to see
if there are any extremely unusual points,
as indicated by the value in the “Contrib.”
(a) These can be selected using the criterion
window, as in step 7 above, but in this
case select “Contrib. >2.0,” etc. to find
(i) To view individual points, select
the two binocular icons to have
the viewers zoom to the selected
tie points to examine particularly
unusual points. If a visual
examination shows that a given
tie-point pair is wrong – i.e. that
the points correspond to entirely
different points on the landscape,
then right-click on the row of the
offending point and select “Delete
selection” to remove it. This will
delete all the points currently
selected, so make sure that only the
desired point is selected.
SOP 2. Preprocessing Landsat Imagery 1
(ii) After deleting a point, recompute
the solution by hitting the sigma
function button again. This
will compute the overall error
without the offending point, and
will redistribute error among
the remaining points. Check for
unusual points again.
(1) This iterative process is done
automatically in the ITPFind
program, so this process likely
should be unnecessary.
(c) Once the reprojection solution has been
computed from the tie points, the input
image can be resampled.
(i) First, record the Control Point
Errors computed in the GCP Tool
window, both for X and Y and for
the total. These should be attached
to the metadata for the resampled
input image. See SOP 5 Data
(ii) Click on the resampled icon in the
“Geo Correction Tools” window
(see fig. 8 above).
(iii) Follow the instructions in section
IV, step 4.b. iii.(3)–(9) above for
(d) When resampling has concluded,
close all windows and do not save the
projection parameters or model. These
would become attached to the images,
which causes confusion if later images
are to be registered to these same
At this point, the input image should be in the same
projection as the geometric-reference image. To check this,
follow these steps.
5. Checking for geometric matching of images.
a. Open a viewer.
b. Open the reference image using the standard file
opening procedure (section IV, steps 1–2 above).
c. In the same viewer, use the open folder icon to start
the “Select Layer to Add” window (fig. 18) and
navigate to the input image and select it, but do NOT
click on OK.
i. Click on the “Raster Options” tab in the Select
Layer to Add window, and uncheck the box next
to “Clear display.”
ii. Click on OK button.
d. Both images should be displayed in the same viewer.
e. In the viewer, select “Utility/Flicker.”
i. The Viewer Flicker dialog will open.
ii. Hit the “Manual filter” button to toggle between
the two images in the viewer.
iii. Alternatively, click the box “Auto Mode” to have
the toggling occur automatically.
iv. When the two images are toggled, there should
be minimal apparent “jumping” of pixels back
and forth across dates. The flickering image
should appear geometrically stable, although
radiometric changes may cause the illusion of
some motion as the images flicker.
v. Scroll around the image and ensure that the
image is stable across the entire range.
vi. Zoom in on some locations and view the flicker
vii. Pay particular attention to mountainous areas.
If the automatic flicker mode is on, if there is a
problem with the terrain correction, the flicker
will reveal distortions that are correlated with
the position of hills and mountains, creating a
pseudo-three-dimensional effect as flickering
(1) If this is the case and terrain correction has
been done in Imagine, this suggests that
either the terrain correction was inaccurate,
or that the original reference image was in
fact NOT terrain corrected.
(a) Check the elevation units of the
DEM and make sure these are in the
appropriate units. If not, redo the terrain
correction of the input image.
(b) Check the history of the geometric
reference image to confirm that it was
terrain corrected or orthorectified.
(c) Repeat the reprojection of the input
image process to the point where tie
points have been loaded in the GCP
tool, and check carefully for unusual tie
points. Check each tie point with high
RMS error and ensure that all points
make sense visually.
Once the input image has been reprojected to the
geometric properties of the geometric reference, and once it
has been confirmed to be accurate, the input image must be
clipped to the study area. Follow the steps in section V to
clip the input image. Use the naming conventions in SOP 5
(Data Management) to name the clipped, terrain-corrected
2 Protocol for Landsat-Based Monitoring of Landscaped Dynamics at North Coast and Cascades Network Parks
VII. Compensating for Sensor and
Atmospheric Influences in Reference
and Input Images
The atmosphere affects the reflectance signal that
impinges on the sensor from the surface of the Earth.
Differences in atmospheric conditions between dates of
imagery will cause differential artifacts in the images that
will confuse later change detection. Therefore, all reasonable
efforts must be made to compensate for these effects. The
first step is to bring both images into a common system of
The units of a Landsat image are simply digital numbers
(DNs), corresponding to the magnitude of energy being
measured by each of the sensor elements in the satellite.
Because the engineering properties of the sensor are known,
the DNs can be converted into physical units of radiance.
Knowing the emission spectrum of the sun entering the
atmosphere, these units of radiance can further be quantified
as a proportion of the incoming radiation ranging from 0 to
1.0. This is known as “top-of-atmosphere reflectance.” It does
not take into account the effects of the atmosphere.
Without direct measurements of atmospheric absorption
at the moment of image acquisition, it is impossible to know
this atmospheric effect directly, but it can be approximated. To
approximate the effects of the atmosphere, it must be assumed
that there are some objects on the surface with little to no
reflectance. These so-called “dark objects” should indicate
a reflectance of zero or near-zero. Because the atmosphere
introduces scattering between the objects and the sensor,
the apparent reflectance of these objects from the sensor’s
perspective is non-zero and positive. The additive effect of
this scattering can be removed by simply subtracting the offset
above these dark objects from all pixels in the image. The
multiplicative effects of the atmosphere can be approximated
simply by a correction factor that scales with the path length
through the atmosphere, which requires only that the sun angle
and elevation be known. All of these steps can be incorporated
in a single-transformation process and placed in a graphic
model in Imagine. Because these were best described in
(Chavez, 1996) as COS-theta or COST methods, these steps
are referred to here as COST processing.
1. Selecting dark-object values for each band:
a. The first part of the COST processing is selection of
dark-object values for each of the six visible bands in
b. An Excel spreadsheet is in the appendixes and serves
as a template for which the relevant values for COST
processing can be stored.
i. This spreadsheet file is named “COST_
ii. Save this Excel file under a new name that
corresponds to the name of the image, so that it
can be connected easily with that image in the
iii. Make a note of this filename in the metadata
associated with this image.
(1) Follow the conventions in SOP 5 Data
Management for naming convention.
c. It is helpful to use the entire image area, not just the
area within the study area of the park, for locating
potential dark objects because the population of
potential targets is much higher. Therefore, for a
given image, go back to the first version of the image
that includes the entire footprint of the original
Landsat image directly after importing into Imagine
i. This image should be the one that was imported
in section III above. It need not be geometrically
referenced for selection of dark-object values.
d. Open the image in a viewer using the Pseudocolor
i. Use the open-folder icon to bring up the “Select
layer to Add” window (shown in fig. 18 above).
(1) Navigate to the desired image folder, select
the image from the available files list so that
its name appears in the “File name:” field.
Do NOT click on OK yet.
ii. Click on the “Raster options” tab.
(1) In the field labeled “Display as:” select
(2) Choose the desired layer. If just starting,
select 1 for band 1, which corresponds to the
blue wavelengths in Landsat.
(3) Retain defaults for the other options.
iii. Click on OK button.
iv. Right-click in the viewer and select “Fit image to
v. Select “Raster/Attributes.”
(1) The Raster Attribute Editor will popup.
(a) There are three columns: Histogram,
Color, and Opacity. The important
columns for this exercise are the
Histogram and the Color columns.
(b) The Histogram column indicates the
count of pixels that have the value
indicated in the row value on the left-
hand side of the attribute editor. For
Landsat images, these range from 0 to
SOP 2. Preprocessing Landsat Imagery
(c) There is a high histogram count in the
row with value 0. These are the pixels
outside the edge of the image area, and
should be ignored.
(d) For bands 1–3 and often 4, the row
number with the lowest non-zero
histogram value is not the row 1, but
rather another row with a higher value.
This indicates that even the darkest
objects in the scene have non-zero
reflectance, and indicates the level of
atmospheric scattering. It is not wise
to simply take the row number of the
lowest non-zero histogram count as the
dark-object value, however, because
these low values could be artifacts
of processing upstream of image
acquisition or even of steps conducted
(i) It is most robust to use images that
only have been resampled using the
NN method, which is why SOP 1
recommends ordering imagery
using the NN option from EDC.
2. Start with the lowest non-zero histogram value and
examine the location of the pixels corresponding to that
value using step 1.a.i.–vii. below, and step sequentially to
higher values until a pixel value is identified whose pixels
reside in the landscape, and which correspond to features
that are expected to be dark-objects step 1.b.i.–ii. below.
a. To examine where pixels at a given pixel value reside
on the image, do the following:
i. Select the row with the desired pixel value. For
example, select the row with the first non-zero
ii. Then place the mouse over the box in the column
labeled “Color” for this row, and right-click.
(1) A small window with a variety of color
options will pop up. Select a color that
allows easy contrast with the grey tones of
the image, perhaps red or yellow.
iii. Look in the viewer at the image and see where
the red pixels reside.
iv. Evaluate these pixels according to the criteria in
step b. below.
v. If these pixels are not considered valid dark
objects, select “Edit/Undo Last Edit” in the
Raster Attribute Editor to change the color patch
back to its original nonhighlighted shade. Then
move to a pixel value one step higher and repeat
vi. If these pixels are considered legitimate dark
objects, then take this value as the approximate
(1) Because few objects are truly nonreflective,
it is typical to assume that the observed
dark objects have an actual reflectance on
the ground of roughly 1 percent. Therefore,
the true value for zero reflectance will be
somewhat lower than the observed value
in the image. The correction factor can be
approximate, since these dark objects are
an approximation of the actual atmospheric
effects. Therefore, take the tentative dark
object value and adjust it as follows to
calculate the actual dark-object value for the
(a) For bands 1–3 and band 7, subtract 1
from the observed value.
(b) For bands 4 and 5, subtract 2 from the
(i) In both cases, however, the
minimum allowed dark-object
value is zero.
vii. Enter this value in the dark-object value column
for the row corresponding to the band in the
Excel spreadsheet saved in step 1.b. above.
b. Evaluating dark-object pixels.
i. Characteristics of valid dark-object-type pixels.
(1) The pixels reside on the landscape, not on
the margins of the scene area.
(2) The pixels are in areas where dark
reflectance is expected:
(a) Water bodies.
(b) Topographic shadow.
(c) Deep shadow on shaded aspects of
(ii) Characteristics of dark objects that are artifacts:
(1) The pixels reside on the very margin of the
active area of the image.
(a) These typically are formed because of
bilinear or cubic convolution of zero
values outside the image with non-zero
values inside the image.
(2) The pixels reside directly next to a very
bright pixel and the image has been
subjected to cubic convolution resampling
Protocol for Landsat-Based Monitoring of Landscaped Dynamics at North Coast and Cascades Network Parks
Documents you may be interested
Documents you may be interested