1. Radiometric Operations

A Radiometric operation uses or changes the values the pixes in an image. In this section, you will apply some of the most common radiometric operations in image processing. But first, you need to get familiar with how sensors capture Electromagnetic radiation and how we can Display an image based on the Tri-stimuli theory. And how we represent colour on paper using the principles of the Subtractive colour system.

1.1. Image Display

Important

Resources. You will require the latest LTR version of LTR version of QGIS, plus the dataset data_histogram_operations.zip. When you unzip the dataset, you will find the following files inside:

• SPOT270611.img – SPOT 5 image with 4 spectral bands.

• tm24aug99.img – TM image of the 24th of august 1999.

• tm24aug99_sub.tif – a subset of tm24aug99.img

• topo34f.img – scan of a topographic map in RGB.

We can enhance the visualisation of images through well-known plugins and providers available in QGIS. In this exercise, we will make use of some of the QGIS tools. Note that none of these tools changes the actual values stored in the raster datasets. They simply change the way the image is displayed to highlight features that are not obvious when using the default visualisation settings.

Task 1

Disable the default contrast stretch. Go to Settings > Options > Rendering tab and scroll down to Contrast enhancement.

Set the algorithms Single band grey, Multiband colour (byte/band) and Multiband colour (> byte/band) to No Stretch. Set the Cumulative pixel count cut to \(2.0\) and \(98.0\%\). Click OK to confirm. See Figure Fig. 1.5

Fig. 1.5 Changing the default settings of contrast enhancement to ‘No Stretch’

1.1.1. Single Band Display & Relative Brightness

In general, remote sensing images can be displayed using Pseudo colour and Colour composite. A common colour composite is the so-called True Colour. Single-band display uses pseudo colour.

Task 2

Create a new QGIS project and open the topo34f.img. If required, change the colour composite for this layer such that Hydrographic elements like water bodies display in cyan colours (use the legend in the right down corner as reference). Right-click on the layer > Properties > Symbology; change the band selection for Red, Green and Blue to Layer_1, Layer_2, and Layer_3.

Task 3

Add the Spot270611.img to the project. Change the display from the default Multiband colour to Singleband Gray and select ‘Band 1’ as the gray band. Right-click on the layer > Properties > Symbology > Render type > Singleband gray > Apply. The band will be displayed in a greyscale with poor contrast. See Figure Fig. 1.6

Fig. 1.6 Displaying band 1 of ‘Spot270611.img’ as greyscale

Then, set min and max values for the contrast stretch. Set contrast enhancement to Stecht to MinMax. Select Cumulative pixel count cut and set the limits to \(35\%\) and \(98\%\). Set Accuracy to Actual (slow). See Figure Fig. 1.7. Click Apply. This will copy the DN values associated with 35 and 98 cumulative percentages the Min and Max of the contrast, respectively.

Fig. 1.7 Contrast enhancement with ‘Stretch to MinMax’

We chose a \(35\%\) for the minimum because the raster file does not contain image data values for3 the whole scene. Approximately \(35%\) of the image includes pixels with a DN Value of 0, which in this case means No Data. After applying the settings above, the image will look like this:

Repeat the previous task. This time apply a MinMax Stretch to all bands of Spot270611.img. You can copy a layer by doing right-click on the ’Spot270611.img’ and then choose Duplicate. Rename each layer name such that it includes the band number, see the example below.

Task 4

Compare the results of each band by toggling the visibility of the layers off and on. Give special attention to the comparison of band 3 and band 2. These two bands are displayed with similar composition, and jet they look quite different from the others. This proves that those bands captured different spectral properties.

Task 5

Use the topographic map topo34f.img to find areas in the Spot270611.img with Water (Cyan), Buildings (Purple) and Evergreen Forest (Green with overprinted symbols). Then, identify the relative brightness in each of the four bands associated with the areas listed above. Fill in the table below.

Do not spend too much time in identifying representative objects and filling the table. Remember that the decision, whether something is grey or light grey is subjective; thus use the same subjectivity when you fill in the table. If you think a class is represented with more than one brightness, you may select more boxes.

Note

Reflection. While working on the previous task, you should have noticed that different classes of land cover may have similar or different brightness within a specific spectral band. Moreover, the same class of land cover may have different brightness in other bands.

1.1.2. Multiband Display: Colour Composites

This section will help you to understand the relationship between the spectral property of a class, the selection of spectral bands for visualisation, and the choice of spectral bands in a colour composite. Suppose you have a SPOT XS image which includes land cover the types: soil, vegetation and water. Such an image will be displayed with a contrast stretch with the band combination of 3, 4 and 2 for RGB.

Attention

Question. Which colours will the land cover types above will have in the colour composite mentioned above? You were introduced in the lecture on how this works for one and two spectral bands. Now it is your turn to predict the case for a 3-bands colour composite.

Use the reflectance curves below to estimate how much each of the three land cover types will reflect relative to each other, for each of the ’SPOT XS’ bands.

For this, assume that there are no other land cover types in the image; i.e. the land cover with the highest reflectance will have the highest Digital Number in a band. Thus, it will be displayed with the highest brightness (i.e., math:DN=255). Likewise, the land cover with the lowest reflectance in a band will be displayed with the lowest brightness (i.e., math:DN=0). For the land cover types that fall between reflectance values, use linear interpolation to estimate its brightness value.

Task 6

Copy RGB brightness values that you estimated into the table. Then use an RGB calculator to determine the approximated colour of each land cover type in this colour composite.

Major land cover	Band 3 (Red)	Band 4 (Green)	Band 2 (Blue)	Colour in composite
Dry-bare soil
Vegetation
Water

Note

Reflection. Based on only theoretical knowledge, you can predict the colour of a land cover type in a colour composite. This is assuming there are no other land cover types present in an image, which may have a lower or higher reflectance. Take the time to understand this statement and study the topic once more if necessary. You could also discuss this topic further in the virtual classroom.

In the previous task, you determined the relative brightness for three types of land cover. Then, you determined the theoretical colour in a specific colour composite. In the following task, you will compare the theoretical colour for water and vegetation with the actual colour in the image when using a colour composite.

Task 7

Open the Spot270611.img in QGIS and use a band combination of 3, 4 and 2 for Red, Green and Blue. Set the contrast enhancement to \(35\%\) and \(98\%\) for all bands. Use Actual (slower) for Accuracy, as shown below.

Zoom in to the image so that you can easily confirm what colour was assigned to water. Also verify whether the theoretical colours, determined in the previous section, match the colours that you see in this colour composite.

Attention

Question. When comparing the theoretical colours and the result of the colour composite, you will find out that there are discrepancies. What could be the courses behind such discrepancies?

Note

Reflection.

You should realise that by knowing the relative spectral reflectance of a class in every spectral band helps to understand and interpret images. Such spectral reflectance is visualised as relative brightness. Moreover, you could predict, for a given sensor, the ranges of DN values for a specific object, if you know and understand the spectral properties of such objects.

Some type of land cover might have similar or the same spectral reflectance property in specific ranges of the Electromagnetic (EM) spectrum. The integration of more bands in the analysis, and the assessment of the spectral properties in other ranges of the EM spectrum, can result in the successful discrimination of more types of land cover. This is true not only in the case of land cover but also for many other objects in a remote sensing image.

---

1.2. Image Enhancement by Histogram Operations

Image enhancement describes a set of operations that aim to improve the way certain features in an image are displayed. One method to enhance images, for interpretation or analysis, uses Histogram operations. These operations use the Histogram of an image to control how it is displayed, and they are usually known as global contrast enhancement operations.

Important

Resources. You will require the same dataset as in the previous section: data_histogram_operations.zip.

• To experiment with contrast enhancement, we will use a TM image of the 24th of august 1999, tm24aug99.img. This image covers a large water body and different types of land cover, and it has areas covered by clouds.

1.2.1. Contrast Stretching

Task 8

Display the tm24aug99.img using the band combination 4, 5 and 3. Set the Stretch to MinMax to a Mean +/- standard deviation of \(2.0\), and the Accuracy to Actual (slower); as shown below.

Task 9

Zoom in to an area covered by clouds over the mainland (centre right). Go to Layer Properties > Symbology. Change the Statistics extent to Current canvas and click OK.

Click on Zoom Full to zoom out to see the whole image. Instead of seeing clouds in whites and pinks, you now can see clouds in several colours.

Fig. 1.8 Top: Image ‘tm24aug99’ with global contrast enhancement. Bottom: Image ‘tm24aug99’ with local contrast enhancement for clouds

Attention

Question. Do you also think that the mainland is displayed in a not-so-good way?

Task 10

Zoom into an area with mostly land and some water, and re-apply the same Contrast enhancement method as before. You can use the context-sensitive option Stretch using current extent as shown below, but verify that it does what you intend.

You will notice a change in contrast in the image. This is because the part of the image that is currently visible include different types of land cover; thus, different local statistics. The results of a contrast stretch based on local statistics changes when the range of values used in the computation changes.

Task 11

Set the Contrast enhancement back to Mean +/- standard deviations and the Statistics extent to Whole raster. Then, right-click on the ’ tm24aug99’ layer and do Export > Save as… For Output mode select Rendered image and enter a self-explanatory name for the output file. Save the file to an existing folder. The file will automatically be added to the Map View.

Attention

Question. Is the exported image different from the original tm24aug99.img?

Compare the properties of the two files, especially the histograms. Determine whether they only look-alike, or if they are the same. You can use the Value tool for the comparison by setting it up in such a way that the tool shows the values of the same band for both images, as shown below.

1.2.2. Choosing Min and Max values

To choose the min and max values for a contrast stretch, the user has to consider which areas of an image are of interest, or which types of land cover are relevant for specific purposes.

Task 12

Remove the exported image from the project; keep only the original image. Zoom into an area on the mainland which is primarily dark orange/brown; they represent areas with forest. Then, on Layers Panel, right-click the original image > Properties > Symbology > Min / Max Value Settings. Set the Statistic extent to Canvas extent. See Figure Fig. 1.9. Click Apply. This will compute the histogram statistics only of the par to the image that it is visible in the Map View. Notice that there are other options to compute statistics than just Min/Max. For example ‘Mean +/- standard deviation’.

Fig. 1.9 Applying local statistics for histogram stretching based on the ‘Canvas extent’.

Note

Reflection. To correctly apply contrast enhancement for specific types of land covers, you need to know which are the types of interest. What their spectral signatures are; the specifications of the spectral bands of the sensor which you have chosen; and you need knowledge of additive colour mixture.

Task 13

Add the tm25aug99_sub.img to the project, and display it using a band combination 4, 5 and 3. this image covers shallow water and land with various types of land cover. Analyse the histograms of the three bands for this image; right-click on the layer > Properties > Histogram. Compute the histogram if necessary. Then, Prefs/Actions > Show selected band and choose the band you want to inspect.

Attention

Question. In which band on display do you expect a significant difference in DN Values between water and land? Use your knowledge on EM radiation. A bi or tri-modal trend in the histogram is an important clue.

Task 14

Use the histograms to identify approximate values for a contrast stretch which will enhance the contrasts between types of coverage on the part of tm25aug99_sub.img with land. Save the result using Export As.. > Rendered image. Remove the resulting layer from the project.

Task 15

Repeat the previous task. This time use the histograms to set a contrast that will enhance the image for distinguishing shallow water.

Note

Reflection. It should be clear to you that for some contrast enhancement methods, statistics of the data play an important role, e.g. mean and standard deviation, minimum and maximum. However, when you know the spectral properties of the objects of interest, the characteristic of the scene, and the sensor; you can interpret the histograms directly and make improvements to make effective use of the brightness values in an image.

---

1.3. Image Enhancement by Filter Operations

Filtering describes a set of radiometric operations used to enhance images. Filters are applied to images for the sake of Noise reduction, Edge Detection, and Edge enhancement.

Important

Resources.

You will require the latest LTR version of LTR version of QGIS, plus the dataset Data_Filter_Operations. When you unzip the dataset, you will find the following files inside:

• tm1999_b4.tif – A scene from band 4 of the Landsat TM. Enschede in 1999.

• tm1999_xs_ml_classification.tif – Classification of the types of land cover in Enschede. From Landsat TM 1999.

• SW-NE_3x3.txt – Definition of a custom filter.

QGIS offers the possibility to apply all kind of filter kernels on images which are displayed in a viewer. In this exercise, we will use tools that apply filters and store the output as temporary files. In such a way, we can easily compare different results.

Task 16

Install the Profile tool plugin. Go to Plugins > Manage Install Plugins, and install the plugin.

Task 17

Set the default contrast stretch to use the \(2 \%\) and \(98 \%\) of the cumulative pixel count for grayscale images. In the Settings menu, go to Options > Rendering tab. Scroll down to Contrast enhancement settings, and set the default value for Single band gray to Stretch To MinMax. Then, set Limits (minimum/maximum) to Cumulative pixel count cut. Make sure that the cut limits are set to \(2.0 \%\) and \(98.0 \%\). Click OK.

Note

QGIS. QGIS is not specifically tailored for Remote Sensing and does not provide filter tools directly. Such filtering tools are available through the processing toolbox using external providers like SAGA and GRASS.

For an overview on how to use the Processing Tools in QGIS, watch the introduction to processing video tutorial.

Task 18

In the Settings menu, go to Options > Processing and check that you have the SAGA and GRASS providers enabled.

1.3.1. Linear Filters

1.3.1.1. Smoothing Filter

Task 19

Apply a linear filter to the ’tm1999_b4’ image. In QGIS, open the tm1999_b4.tif. Your project should assume the same Spatial Reference System as the image (EPSG:32632 WGS84/UTM zone 32N). In the Processing Toolbox, open the SAGA tool called User defined filter.

Confirm that ’tm1999_b4’ is the input and click Default filter matrix (3x3) to open an empty filter kernel. Enter the weights of an average filter kernel. Ensure that the sum of weights is equal to 1. Confirm with OK. In the User defined filter dialogue execute the kernel by clicking OK. The output is added to the Map View as a temporary file.

In the Layers panel right-click the ’Filtered Grid’ layer and rename to ’Average’.

Attention

Question. Which kernel weights did you use in the previous task? Write them down.

Task 20

Explore the filter results around the Twente Airport. Reset the zoom to fit the image to the Map View. Next, change the scale, in the box at the bottom of the Map View, to \(1:75,000\).

Zoom in to the major runway of the Twente Airport. See Figure Fig. 1.10.

We will use the Profile tool to compare the results of the average filter and original image. If you do not know how to install the Profile Tool plugin, watch the video tutorial on installing plugins in QGIS.

Start the Profile tool. Select one of the layers in the Layer Panel. Zoom into the centre of the image and draw a profile (line) across the major runway. Click on Add Layer. Select the other layer in the Layers panel and click on Add Layer again.

Fig. 1.10 Comparison of the result of a smoothing filter to ‘tm1999_b4.img’ using the ‘Profile tool’

Attention

Question. Just by looking at the graph of the Profile tool, can you tell which profile belongs to the layer with the average filter?

Task 21

Draw profiles at different locations, and confirm your knowledge of the effects of applying an average filter (smoothing) to the image.

1.3.1.2. Gradient Filter

Task 22

Use the User defined filter tool to apply a filter using the weights in the figure below, to the original ‘tm1999_b4’ layer. Rename the resulting layer to ‘Laplace’.

Error

You will notice that the output of the Laplace filter contains many pixes with No Data (hollow pixels). Possibly, you also noticed an error message when executing the filter. This is an issue with the tool itself, possibly a bug. Please, ignore the No Data pixels this time and concentrate on the outcome in the rest of the image.

Attention

Question.

• Is the kernel above a detection kernel? If yes, what does it detect?

• Does the layer resulting from the previous task contain the same brightness as the original image for area objects?

Task 23

Examine the result of the Laplace filter. Toggle on and off the visibility of the ‘Laplace’ layer to check what has happened. Zoom in to the edge of the image until you see individual pixels. Toggle the visibility of the ‘Laplace’ layer again.

Attention

Question. What phenomena do you observe? Can you explain it?

Task 24

Open the histogram of the ‘Laplace’ layer. Go to Properties > Histogram; check the values in the image.

Attention

Question. Around which value does the histogram has its centre?

Task 25

In the Profile tool add the ’Laplace’ layer and toggle the visibility for the other layers. Confirm that the filter kernel detected two edges, i.e. both sides of the runway.

The ’Laplace’ layer looks rather artificial. The brightness of the original image is gone; the lighter and darker areas in the original have now a common grey tone and high contrasting pixels at the edges. This filter has detected the changes (edges) between local lighter and darker pixels. The circular build of the kernel that you applied, i.e. all negative weights around the centre with positive weight, detected changes in all directions.

Note

Reflection. We can use the output values of an edge detection kernel to discriminate edges with high contrast and edges with low contrast. For example, a field with relatively low values and a neighbouring field with relatively high values will result in edges with high contrast. In comparison, edges with low contrast may be the result of adjacent pixels which have a slightly different value. You could use a threshold to select only edges with high contrast and delimit the edges between areas or fields.

Directional filter kernels have positive and negative weights, and their output values could be positive or negative. Positive values represent edges which correlate positively with the kernel, negative values represent edges which correlate negatively – ‘ they are opposite to’ –. We can also use the sign of the values to discriminate different classes or types of edges.

If you do not understand these statements, review the lecture material once more.

1.3.1.3. Edge Enhancement

Task 26

Repeat the steps of the previous task, but this time use kernels with the following values for the centre: 12, 16 and 200. Increasing the centre value will increase the weight of the centre pixel in the original image. When using a value of about 16 for the centre pixel; the kernel will calculate the Laplace enhancement of the image. Then, the resulting layer will look like the original image.

1.3.1.4. Custom Filters

The Laplace kernel detects edges in all directions. We can also define kernels which detect edges in specific directions. In the next task, you will use a costume filter defined in a text file. Check the content of SW-NE_3x3.txt by opening in a text editor.

Attention

Question. Which are the weight of the filter define in SW-NE_3x3.txt? Write them down.

Task 27

In the Processing Toolbox, use the r.mfilter tool of GRASS to detect edges in a specific direction. Use the file SW-NE_3x3.txt as Filter file.

Check the results and confirm that one of the runways of the Twente Airport was not detected at all! Also, ensure that the edges of the main runway were detected. The result should show positive values on one side of the runway and a negative on the other. This is because of the correlation of the results with the positive and negative weights in the kernel.

Notice also, that the edges in the resulting layer have a slope and a direction. In this case, the kernel has detected edges in the SW-NE direction.

1.3.2. Enhancement using Non-linear Filters

1.3.2.1. Rank-Order Filter

Task 28

In the Processing Tools, look for the SAGA Rank filter. Select the ’tm1999_b4’ layer as the input grid. Set the Search Mode to Square and a Radius of 1. Use a Rank (Percent) of 50. This settings essentially define a median filter. Execute the filter.

Use the Profile tool to inspect the results and confirm the difference between the original image and the result of the Average and Median filters. Pay special attention to locations where you expect variations, for instance, around the edges between areas.

1.3.2.2. Majority Filter

Task 29

Add the tm1999_tm_xs_ml_classification.tif into the Map View; this layer contains a set of classes representing land cover. You will notice that there are many isolated pixels inside some homogeneous areas. For example, the yellow pixels identify maise, which usually does not grow in such small parcels.

In the Processing Toolbox, open the r.neighbors tool. Select ’tm1999_xs_ml_classifciation’ as input raster and set the neighbourhood operation as ’mode’ (also known as majority filter).

The results will be shown as greyscale or black and white. You can apply a pseudocolour to visualise the results properly. Copy the style of the ’tm1999_xs_ml_classifciation’ into the result of the majority filter layer. See Fig. 1.11

In the Layer panel, right-click on the ’tm1999_xs_ml_classifciation’ layer and copy the style. Then, right-click on the majority filter layer and paste the style.

Fig. 1.11 Copying the style between two raster layers in QGIS

Task 30

Toggle the visibility of the filtered result and compare it with the original landcover layer. Confirm that most isolated pixels have disappeared, and that thin lines of pixels surrounded by homogeneous areas also disappeared (e.g. a runway in the Twente Airport). The main runway should still be distinguishable.

Attention

Question. Can you explain why one of the small runways of the Twente Airport disappeared after applying a majority filter?

Task 31

Use the r.neighbor to compute another mode filter on the ’tm1999_xs_ml_classifciation’ layer. This time use a neighbourhood of size \(5\); which is to apply a \(5 \times 5\) kernel. Compare the result of this filter with the \(3 \times 3\) majority filter.

Attention

Question. What do you observe when comparing results of \(3 \times 3\) and \(5 \times 5\) majority filters around the main runway of the Twente Airport?

Task 32

Experiment with the application of consecutive filters. Apply a \(3 ~\times~ 3\) majority filter to the results of the already existing ‘3x3 majority’ filter layer. Compare these results with the result of applying a single \(5 \times 5\) majority filter. You should see that the results are not the same.

Note

Reflection. In summary. You should acknowledge that in the case of neighbour operations, such as the majority filter. The results change depending on the size of the kernel, and the number of times a filter is applied to a dataset.

Section author: Wan Bakx, André Mano & Manuel Garcia Alvarez