# PipelineCropAndMaskFromRegion.sml

# 11 February 2008
# Randy Smith, MicroImages, Inc.

# Sample standalone script to demonstrate use of a region as a pipeline
# source to mask and crop a source image.

# Script inputs: a single raster in a Project File, and a vector object
# with polygons designating the part of the input raster to retain.
# The polygons are converted to a region, which is then used to create a
# pipeline region source.

# Script output: a new raster with areas outside the polygons masked out,
# and cropped to the extents of the polygons.

# Requires Version 2008:74 or later of the TNT products.

numeric err;

# error checking procedure
proc ReportError(numeric linenum, numeric err) {
	printf("FAILED -line: %d, error: %d\n", linenum - 1, err);
	PopupError(err); 
	}


################################################################################
#################################  Main program  ###############################
clear();

# set context so script does not exit on error, allowing manual
# handling of errors.
_context.AbortOnError = 0;


# CHOOSE INPUT RASTER
class RVC_OBJITEM rastInObjItem;
DlgGetObject("Choose input raster:", "Raster", rastInObjItem, "ExistingOnly");


# SET UP PIPELINE SOURCE FOR IMAGE
class IMAGE_PIPELINE_SOURCE_RVC sourceImg(rastInObjItem);
err = sourceImg.Initialize();
if (err < 0 )
	ReportError(_context.CurrentLineNum, err);
else {
	numeric numLins, numCols;
	numLins = sourceImg.GetTotalRows();
	numCols = sourceImg.GetTotalColumns();
	print("Source image initialized.");
	printf("Image size = %d lines, %d columns\n", numLins, numCols);
	}

# Get georeference and coordinate reference system from the image.
class IMAGE_PIPELINE_GEOREFERENCE georefImg;
georefImg = sourceImg.GetGeoreference();

class SR_COORDREFSYS imgCRS = georefImg.GetCRS();

if (!imgCRS.IsDefined() ) {
	print("Image is not georeferenced; exiting now.");
	Exit();
	}
else if (imgCRS.IsLocal() ) {
	print("Image has local engineering georeference; exiting now.");
	Exit();
	}
else
	printf("Image coordinate reference system = %s\n", imgCRS.Name );

# Get coordinate transformation from image to its defined CRS
class TRANS2D_MAPGEN imgToCRS = georefImg.GetTransGen();


# CHOOSE VECTOR OBJECT WITH POLYGONS TO INDICATE AREAS TO PROCESS
class RVC_VECTOR polyVect;
class RVC_OBJITEM vectObjItem;
class RVC_GEOREFERENCE vectGeoref;
class SR_COORDREFSYS vectCRS;
DlgGetObject("Choose vector object with polygons outlining areas to process:", "Vector", vectObjItem, "ExistingOnly");
polyVect.Open(vectObjItem, "Read");

if (polyVect.$Info.NumPolys == 0) {
	print("Vector object contains no polygons; exiting now.");
	Exit();
	}

err = polyVect.GetDefaultGeoref(vectGeoref);
if (err < 0 ) {
	ReportError(_context.CurrentLineNum, err);
	print("Polygon vector must be georeferenced and is not; exiting now.");
	Exit();
	}

vectCRS = vectGeoref.GetCoordRefSys();
printf("Vector coordinate reference system = %s\n", vectCRS.Name);


# CHOOSE OUTPUT RASTER
# DlgGetObject creates only an ObjItem for the output; the new raster will then be created
# by the TARGET_RVC with the desired compression.  The ObjectPath of the ObjItem is not
# valid until the output raster is actually created by the pipeline.
class RVC_OBJITEM rastOutObjItem;
DlgGetObject("Choose new raster for result:", "Raster", rastOutObjItem, "NewOnly");


# CONVERT INPUT VECTOR TO A REGION FOR USE AS A PIPELINE SOURCE
class REGION2D RegFromVect = ConvertVectToRegion(polyVect, GetLastUsedGeorefObject(polyVect) );

# if vector has different coordinate reference system from image, reproject region
if (vectCRS.Name <> imgCRS.Name) {
	printf("Converting region to image CRS: %s\n", imgCRS.Name);
	RegFromVect.ConvertTo(imgCRS);
	}

class RECT vecRegionExtents = RegFromVect.Extents;
printf("Minimum map extents of polyVect region: %.2f, %.2f\n", vecRegionExtents.x1, vecRegionExtents.y1);
printf("Maximum map extents of polyVect region: %.2f, %.2f\n\n", vecRegionExtents.x2, vecRegionExtents.y2);


# GET EXTENTS OF THE IMAGE SOURCE AS A REGION
class REGION2D imgReg;
err = sourceImg.ComputeGeoreferenceRegion(imgReg);
if (err < 0 )
	ReportError(_context.CurrentLineNum, err);


# CHECK THAT REGION FROM VECTOR IS CONTAINED WITHIN IMAGE REGION
if (imgReg.TestRegion(RegFromVect, "FullInside") ) {
	print("Image region contains vector extents.");
	}
else {
	print("Input vector is not contained within the source image.  Exiting now.");
	Exit();
	}

# PIPELINE SOURCE FOR REGION CREATED FROM THE VECTOR POLYGONS
# Construct using TRANS2D_MAPGEN with coordinate transformation from source image to its CRS.
# Use source image as a reference so the "image dimensions" of the region source match,
# enabling the region source to be used as a mask in the MASKVALIDITY filter
class IMAGE_PIPELINE_SOURCE_REGION sourceReg(RegFromVect, imgToCRS, sourceImg);
err = sourceReg.Initialize();
if (err < 0 )
	ReportError(_context.CurrentLineNum, err);
else {
	print("Initialized region source.");
	printf("Size of region = %d lines, %d columns\n", sourceReg.GetTotalRows(), sourceReg.GetTotalColumns() );
	}

# Get the extents of the region as a RECT in image coordinates
# to use to crop the source image.  These are the extents of the region used to
# create the region source, not those of the image used as its reference.
class RECT regExtents;
regExtents = sourceReg.GetRegionExtents();
printf("Minimum values of cropping rectangle from polygons: x = %.2f, y = %.2f\n", regExtents.x1, regExtents.y1);
printf("Maximum values of cropping rectangle from polygons: x = %.2f, y = %.2f\n", regExtents.x2, regExtents.y2);


# SET UP PIPELINE FILTER TO MASK PORTIONS OF SOURCE IMAGE
# mask area outside the region created from the vector polygons
class IMAGE_PIPELINE_FILTER_MASKVALIDITY filterMask(sourceImg, sourceReg);
err = filterMask.Initialize();
if (err < 0 )
	ReportError(_context.CurrentLineNum, err);
else {
	print("Initialized image mask filter.");
	}


# SET UP PIPELINE FILTER TO CROP THE SOURCE IMAGE
class IMAGE_PIPELINE_FILTER_CROP	filterCrop(filterMask, regExtents);
err = filterCrop.Initialize();
if (err < 0 )
	ReportError(_context.CurrentLineNum, err);
else {
	print("Initialized image crop filter.");
	}


# PIPELINE TARGET
class IMAGE_PIPELINE_TARGET_RVC target(filterCrop, rastOutObjItem);
err = target.Initialize();
if (err < 0 )
	ReportError(_context.CurrentLineNum, err);
else {
	print("Initialized target.");
	}

target.Process();

print("Done.");