PCI Geomatics Help Center

How can we help you today?

Object Analyst - Batch Classification in Python - Geomatica Banff SP2

PCI Geomatics -

In Geomatica Banff standalone Object Analyst algorithms were added in both Python and EASI. These new algorithms allow users to easily run an object-based classification on many datasets. The first steps in the workflow are to segment, calculate attributes and collect/import training sites for the initial image. Once that is complete the attributes fields are exported to a text file, the training model is generated and finally SVM/Random Trees classification is run on that initial image. Each of the additional images is processed in a batch Python script, using the attribute field name text file and training model that were generated from the initial image. The diagrams below shows the basic workflow.

Additional resources for creating Python scripts using Geomatica algorithms are available from the Developer Zone. Note that the script in this tutorial makes use of Python 3.5 as this is the new default version used in Geomatica Banff.

Process Initial Image – Python

Import necessary modules & setup inputs/outputs

OASEG

OACALCATT

OAGTIMPORT

OAFLDNMEXP

OASVMTRAIN

OASVMCLASS

OARTTRAIN/OARTCLASS

Process additional images in batch

OASEG

OACALCATT

OASVMCLASS

OARTCLASS

Edit, Save and Apply Representation

Edit Classification Representation

Save Representation to RST file

Apply Representation to Additional Vector Layers

Optional: Process Initial Image - Focus

Segmentation

Attribute Calculation

Collect Training Sites

Full Script

SVM training and classification workflow

mceclip0.png

Random Trees training and classification workflow

mceclip1.png

The workflow outlined in the tutorial shows how to fully process the initial image and additional images in Python. Note that you must have at least Geomatica Banff SP2 installed to run this full process including the new algorithm OAGTIMPORT. You can optionally run the segmentation (OASEG), attribute calculation (OACALCATT) and training site collection (OAGTIMPORT) steps in Focus and then continue the python process from OAFLDNMEXP. The steps in Focus are outlined in the section Optional: Process Initial Image - Focus.

For this tutorial we will use Landsat-8 data of southern Ontario, acquired in June over four different years. 

Download Tutorial Data: https://pcigeomatics.sharefile.com/d-s34d5345eb8441d6a 

June 3, 2015 – Initial Image

June 21, 2016

LC08_L1TP_018030_20150603_20170226_01_T1

LC08_L1TP_018030_20160621_20180130_01_T1

 mceclip1.png  mceclip2.png

June 11, 2018

June 14, 2019

LC08_L1TP_018030_20180611_20180615_01_T1 

 LC08_L1TP_018030_20190614_20190620_01_T1

 mceclip3.png  mceclip4.png

Images should have similar land cover types in order for the batch classification to produce good results for each image.

Process Initial Image – Python

The first step is to process the initial image. In this case we will use LC08_L1TP_018030_20150603_20170226_01_T1 (June 3, 2015) as the initial image. The following algorithms will be run on the initial image to generate the training model and complete the classification: OASEG (OASEGSAR), OACALCATT (OACALCATTSAR), OAGTIMPORT, OAFLDNMEXP, OASVMTRAIN (OARTTRAIN) and OASVMCLASS (OARTCLASS). You can check the Geomatica Help for additional information on the required parameters for each algorithm.

Import necessary modules & setup inputs/outputs

The first step in the Python script is to import the required modules and set up the variables that will point to our input and output directories and files. You can import the PCI algo module which will be used to access the Geomatica algorithms. Also import the os module to access the operating system’s file structure and glob module to find files in a directory.

Create input and output variables as outlined below.

from pci import algo
import os
import glob

# Input Variables
# Initial Image - This image will be used to generate the training model
init_image = r"C:\Tutorials\OA_Batch\Initial_Image\June2015.pix"
# Intial segmentation containing the training sites
init_seg = r"C:\Tutorials\OA_Batch\Python_Classifications\June2015_seg.pix"
# Ground truth point file
ground_truth = r"C:\Tutorials\OA_Batch\ground_points.pix"
# Additional Images - The batch classification will be run on these images
add_images = r"C:\Tutorials\OA_Batch\Additional_Images"

# Output Variables
# Output file location
output_folder = r"C:\Tutorials\OA_Batch\Python_Classifications"
# Text file containing exported attribute names
fld = r"C:\Tutorials\OA_Batch\Python_Classifications\att_fld.txt"
# Training model
training_model = os.path.join(output_folder, "training_model.xml")

OASEG

OASEG applies a hierarchical region-growing segmentation to image data and writes the resulting objects to a vector layer. The input file (fili) is the initial image variable (init_image) that we defined above. The output image will be the new segmentation file (init_seg). The scale parameter (segscale) will be set to 35 – this will create larger segments then the default of 25.

algo.oaseg(fili=init_image, filo=init_seg, segscale=[35], segshape=[0.5], segcomp=[0.5])

OACALCATT

OACALCATT calculates attributes of objects (polygons) in a vector layer, and then writes the output values to the same or a new segmentation vector layer.

The input image (fili) will be set to the initial image (image) that is being processed in the loop. Only channels 2-6 (B, G, R, NIR, SWIR) will be used in the attribute calculation (dbic). These attributes are calculated for each of the polygons in the segmentation layer created in OASEG (filv, dbvs). The calculated attributes are then saved back to the same segmentation layer (filo, dbov). The Maximum and Mean statistical attributes (statatt) and Rectangularity geometrical attribute (geoatt) will be calculated. Additionally, NDVI vegetation index will be calculated (index) which will help to classify vegetation.

algo.oacalcatt(fili=init_image, dbic=[2, 3, 4, 5, 6], filv=init_seg, dbvs=[2], filo=init_seg, dbov=[2], statatt="MAX, MEAN",
geoatt="REC", index="NDVI")

OAGTIMPORT

Before the release of Geomatica Banff SP2, OAGTIMPORT did not exist and training sites had to be manually collected in Focus. You can still perform the initial training site collection in Focus as outlined in the section Optional: Process Initial Image - Focus. With the addition of OAGTIMPORT you can now import ground truth points as training sites directly in python.

OAGTIMPORT will import a series of training or accuracy sample points into the attribute table of a segmentation file. The training or accuracy samples must be vector points (not polygons) with a field containing the class names. In this example the ground truth dataset includes points that were collected from the initial image in Focus. You may have existing ground truth points for your own datasets. You can also convert existing training site or classification data from polygons to points using the POLY2PNT algorithm.

The point file (ground_truth) will be set as gtfili. In this point file there is a field, gt_class, which lists the classification associated to this point. This field will be set for gtfldnme. Parameters filv,  dbvs, filo and dbov will be set to the initial segmentation (init_seg) file and segment. The samptype is set to “Training” as we want these points to generate training sites, instead of accuracy assessment sites. Resrule is set to “First” so the first-encountered point is used to set the training class of the segmentation polygon.

algo.oagtimport(gtfili=ground_truth, gtfldnme="gt_class", filv=init_seg, dbvs=[2],filo=init_seg, dbov=[2],
samptype="Training", resrule="First")

OAFLDNMEXP

OAFLDNMEXP exports the names of attribute fields from an OA segmentation file to a text file. The input file (filv) is the initial segmentation file that we created in Focus. dbvs will be set to 2 as this is the segment number of the segmentation vector layer. In this tutorial we will export all attribute-field names from Object Analyst (“ALL_OA”) but fldnmflt can also be set to specific fields. The output text file (tfile) will be set to the fld variable that we created earlier in the script.

algo.oafldnmexp(filv=init_seg, dbvs=[2], fldnmflt="ALL_OA", tfile=fld)

OASVMTRAIN

OASVMTRAIN uses a set of training samples and object attributes, stored in a segmentation attribute table, to create a series of hyperplanes, and then writes them to an output file containing a Support Vector Machine (SVM) training model. The same segmentation file that we created in Focus will be used as the input (filv). The field name text file (fld) will be used as another input (tfile). The kernel function for the SVM classification (kernel) will be the default radial-basis function (RBF). Finally the training model (trnmodel) parameter will be set to the training_model variable we created earlier.

algo.oasvmtrain(filv=init_seg, dbvs=[2], trnfld="gt_class", tfile=fld, kernel="RBF", trnmodel=training_model)

OASVMCLASS

OASVMCLASS uses Support Vector Machine (SVM) technology to run a supervised classification based on an SVM training model you specify. This is the final algorithm that needs to be run on the initial image in order to complete the classification. When this algorithm is run, new fields are added to the output file (filo) which include the classification information. The initial segmentation (filv), field name file (tfile) and training model (trnmodel) from the earlier algorithms will be used as inputs. The output vector file (filo) will be set to the int_seg variable and the output vector segment number (dbov) will be 2, as we want the classification fields added to the initial segmentation vector segment. This will ensure that all object analyst fields (attributes, training, and classification) are contained in a single vector layer.

algo.oasvmclass(filv=init_seg, dbvs=[2], tfile=fld, trnmodel=training_model, filo=init_seg, dbov=[2])

OARTTRAIN/OARTCLASS

Similar to the SVM classification algorithms (OASVMTRAIN/OASVMCLASS) you can specify the Random Trees (RT) classification algorithms instead. RT belongs to a class of machine learning algorithms which does ensemble classification. The classifier uses a set of object samples that are stored in a segmentation-attribute table. The attributes are used to develop a predictive model by creating a set of random decision-trees and by taking the majority vote for the predicted class. 

You can choose to run the SVM algorithms or the RT algorithms. Below is code that includes the RT algorithms.

algo.oarttrain(filv=init_seg, dbvs=[2], trnfld="gt_class", tfile=fld, trnmodel=training_model)
algo.oartclass(filv=init_seg, dbvs=[2], trnmodel=training_model, filo=init_seg, dbov=[2])

Process additional images in batch

Now that the initial image is classified, we can apply that same training model to the additional images.

The first step is to create a list of all valid additional images in the add_images folder. The glob module can be used to create a list of all .pix files in the add_images folder.

file_list = glob.glob(os.path.join(add_images, "*.pix"))

You can iterate through the file_list using a for loop. The algorithms within the loop will be run on each individual image.

A new segmentation pix file (add_seg) will need to be created for each input image. The os module is used to establish the naming for the new segmentation file:

os.path.join(output_folder, os.path.basename(os.path.splitext(init_image)[0]) + '_seg.pix')

  • os.path.splitext() splits the file pathname (…\Additional_Images\June2016.pix) at the period of the extension. The first section [0] of the split is kept -…\Additional_Images\June2016
  • os.path.basename() gets only the basename (not path) of the file - June2016.
  • + ‘_seg.pix” adds this suffix to the split basename - June2016_seg.pix
  • path.join() joins the two path components: directory: output folder and basename: June2016_seg.pix to create a full file pathname - …\Python_Classifications\June2016_seg.pix
for image in file_list:
print("Currently processing:", os.path.basename(image))
add_seg = os.path.join(output_folder, os.path.basename(os.path.splitext(image)[0]) + '_seg.pix')

OASEG

OASEG applies a hierarchical region-growing segmentation to image data and writes the resulting objects to a vector layer.

The input image (fili) will be the current image (image) that is being processed in the loop. The output image will be the new segmentation file (add_seg). The scale parameter (segscale) will be set to 35 – this will create larger segments then the default of 25.

    algo.oaseg(fili=image, filo=add_seg, segscale=[35], segshape=[0.5], segcomp=[0.5])

OACALCATT

OACALCATT calculates attributes of objects (polygons) in a vector layer, and then writes the output values to the same or a new segmentation vector layer.

The input image (fili) will be set to the current image (image) that is being processed in the loop. Only channels 2-6 (B, G, R, NIR, SWIR) will be used in the attribute calculation (dbic). The attributes are calculated for each of the polygons in the segmentation layer that we created in OASEG (filv, dbvs). The calculated attributes are then saved back to the same segmentation layer (filo, dbov). The Maximum and Mean statistical attributes (statatt) and Rectangularity geometrical attribute (geoatt) will be calculated. Additionally, NDVI vegetation index will be calculated (index) which will help to classify vegetation.

    algo.oacalcatt(fili=image, dbic=[2,3,4,5,6], filv=add_seg, dbvs=[2], 
filo=add_seg, dbov=[2], statatt="MAX, MEAN", geoatt="REC", index="NDVI")

OASVMCLASS

We will now run OASVMCLASS in a similar method to how we ran it on the initial image. The same training model from OASVMTRAIN and attribute field text file from OAFLDNMEXP will be used. The current segmentation file (add_seg) will be used as both the input and output (filv/dbvs, filo/dbov).

    algo.oasvmclass(filv=add_seg, dbvs=[2], tfile=fld, trnmodel=training_model, filo=add_seg, dbov=[2])

OARTCLASS

If you ran OARTTRAIN/OARTCLASS when processing the original image you will want to use that same algorithm for the additional images.

    algo.oartclass(filv=add_seg, dbvs=[2], trnmodel=training_model, filo=add_seg, db

Once the batch processing is complete each segmentation file will include the classification fields – Label (integer class label), Class (string class label), Prob (class voting probability). You can then open the vector files in Focus and apply a representation to view the classification.ov=[2])

Edit, Save and Apply Representation

When a segmentation vector layer is loaded into Focus you will need to adjust the representation of the vector layer to show the classification. This representation can then be saved to an RST file and that RST file can be used to load the same representation for the additional vector layers.

Edit Classification Representation

The following steps outline how to adjust the representation of the first segmentation vector file in Focus to display the classification results.

  1. Open the first segmentation pix file in Focus
  2. Right-click the segmentation vector layer in the Files tab and choose View
  3. In the Maps tab, right-click the new vector layer > Representation Editor
  4. In the representation editor, change the Attribute option to the supervised classification field
  5. Click More >>
  6. Make sure the method is Unique Values is selected
  7. Make sure the Generate new styles option is selected and that you choose the style that you want to use.
  8. Click Update Styles
  9. You can then change the colour of each class in the top section of the panel.

mceclip0.png

Save Representation to RST file

In order to apply the same representation to various files you will need to first save the representation.

  1. In Focus > Maps tab > right-click the new vector layer > Representation Editor
  2. To save a representation click mceclip1.png on the representation editor.

mceclip2.png

Apply Representation to Additional Vector Layers

When you load the additional classifications to Focus you want to ensure that the RST file you created is linked to the map in your project. This representation can then be easily applied to the additional images.

  1. In Focus > Maps tab > right-click on the Unnamed Map > Representation > Load.
  2. Select the RST file that you just created and click Open.
  3. The RST file will now be listed in the Maps tab

 mceclip3.png

  1. Open the additional segmentation files in Focus.
  2. In the Maps tab, right-click on one of the new vector layers > Representation Editor
  3. In the Representation Editor, change the Attribute drop-down menu to the SupClass field and click OK

 mceclip4.png

  1. The representation from the RST file will be applied.
  2. You can follow steps 5 & 6 for each of the additional segmentation vector layers.

Optional: Process Initial Image - Focus

The workflow above outlines how to run the entire Object Analyst batch process in Python. The example below shows Segmentation, Attribute Calculation and Training Site Collection in the Focus GUI. You can run OASEG and OACALCATT on the initial image before collecting training sites in Focus. 

The full workflow for processing an image in the OA Focus GUI is outlined in the Object Analyst tutorial.

Segmentation

The first step in the Object Analyst workflow is segmentation.

  1. In Focus open the Object Analyst window from the Analysis drop-down menu.
  2. Change the Operation to Segmentation
  3. Click Select beside the Source Channels box.
  4. In this window open the initial image - June2015.pix and select channels 2-6.
  5. Click OK

 mceclip5.png

  1. In the main OA window change the scale to 35 – this will create larger segments then the default of 25
  2. Set an output file and layer name
  3. Click Add and Run

mceclip6.png

Attribute Calculation

The next step in the workflow is to calculate attributes for each of the segments created in the previous step.

  1. Switch the Operation to Attribute Calculation
  2. Keep the same channels checked off in the Source Channels list
  3. In the Attributes to Calculate section select Statistical: Max & Mean, Geometrical: Rectangularity and Vegetation Indices: NDVI

mceclip7.png

Collect Training Sites

For this tutorial ground truth points are supplied which will be imported and used as training sites. More information on collecting training sites is available from the Object Analyst Tutorial.

  1. Change the Operation to Training Site Editing
  2. Make sure that the Segmentation layer is checked off under Training Vector Layer

 mceclip8.png

  1. Click on the Import… button
  2. Select the ground truth file by clicking Browse
  3. In the Layer section select the vector segment GT:Training_Ground_Truth
  4. Select the Training field and Sample type: Training
  5. Click Import

 mceclip9.png

 mceclip10.png

  1. Click on the Edit…button to open up the Training Sites Editing window
  2. Ensure that Training Field is set to Training
  3. You can now view the imported training sites

mceclip11.png

mceclip12.png

Full Script

from pci import algo
import os
import glob

# Input Variables
# Initial Image - This image will be used to generate the training model
init_image = r"C:\Tutorials\OA_Batch\Initial_Image\June2015.pix"
# Intial segmentation containing the training sites
init_seg = r"C:\Tutorials\OA_Batch\Python_Classifications\June2015_seg.pix"
# Ground truth point file
ground_truth = r"C:\Tutorials\OA_Batch\ground_points.pix"
# Additional Images - The batch classification will be run on these images
add_images = r"C:\Tutorials\OA_Batch\Additional_Images"

# Output Variables
# Output file location
output_folder = r"C:\Tutorials\OA_Batch\Python_Classifications"
# Text file containing exported attribute names
fld = r"C:\Tutorials\OA_Batch\Python_Classifications\att_fld.txt"
# Training model
training_model = os.path.join(output_folder, "training_model.xml")

# Process Initial Image
print("Processing initial image:", os.path.basename(init_image))

# OASEG (OASEGSAR) - Segment an image
algo.oaseg(fili=init_image, filo=init_seg, segscale=[35], segshape=[0.5], segcomp=[0.5])

# OACALCATT (OACALCATTSAR) - Calculate object attributes
algo.oacalcatt(fili=init_image, dbic=[2, 3, 4, 5, 6], filv=init_seg, dbvs=[2], filo=init_seg, dbov=[2], statatt="MAX, MEAN",
geoatt="REC", index="NDVI")

# OAGTIMPORT - Import ground-truth points into an existing segmentation
algo.oagtimport(gtfili=ground_truth, gtfldnme="gt_class", filv=init_seg, dbvs=[2],filo=init_seg, dbov=[2],
samptype="Training", resrule="First")

# OAFLDNMEXP - Export names of attribute fields from Focus Object Analyst to a text file
algo.oafldnmexp(filv=init_seg, dbvs=[2], fldnmflt="ALL_OA", tfile=fld)

# SVM Classification
# OASVMTRAIN - Object-based SVM training
algo.oasvmtrain(filv=init_seg, dbvs=[2], trnfld="gt_class", tfile=fld, kernel="RBF", trnmodel=training_model)

# OASVMCLASS - Object-based SVM classifier
algo.oasvmclass(filv=init_seg, dbvs=[2], trnmodel=training_model, filo=init_seg, dbov=[2])

# Random Trees Classification
# OARTTRAIN - Object-based Random Trees training
#algo.oarttrain(filv=init_seg, dbvs=[2], trnfld="gt_class", tfile=fld, trnmodel=training_model)

# OARTCLASS - Object-based Random Trees classifier
#algo.oartclass(filv=init_seg, dbvs=[2], trnmodel=training_model, filo=init_seg, dbov=[2])

# Apply training model and classify additional image in batch
print("Processing additional images in batch...")
file_list = glob.glob(os.path.join(add_images, "*.pix"))

for image in file_list:
print("Currently processing:", os.path.basename(image))
add_seg = os.path.join(output_folder, os.path.basename(os.path.splitext(image)[0]) + '_seg.pix')

# OASEG (OASEGSAR) - Segment an image
algo.oaseg(fili=image, filo=add_seg, segscale=[35], segshape=[0.5], segcomp=[0.5])

# OACALCATT (OACALCATTSAR) - Calculate object attributes
algo.oacalcatt(fili=image, dbic=[2,3,4,5,6], filv=add_seg, dbvs=[2], filo=add_seg, dbov=[2], statatt="MAX, MEAN",
geoatt="REC", index="NDVI")

# OASVMCLASS - Object-based SVM classifier
algo.oasvmclass(filv=add_seg, dbvs=[2], trnmodel=training_model, filo=add_seg, dbov=[2])

# OARTCLASS - Object-based Random Trees classifier
#algo.oartclass(filv=add_seg, dbvs=[2], trnmodel=training_model, filo=add_seg, dbov=[2])
Have more questions? Submit a request

Comments

Powered by Zendesk