# PipelineResampleToMatch.sml

# Sample standalone script to demonstrate a simple image pipeline application:
# resampling/reprojecting a source raster to match the extents, cell size, and
# projection of a reference raster.

# 21 December 2007
# revised 13 February 2008

# Randy Smith, MicroImages, Inc.

# Requires Version 2007:73 or later of the TNT products.

numeric err;														# error code returned

# error checking procedure
proc ReportError(numeric linenum, numeric err) {
	printf("FAILED -line: %d, error: %d\n", linenum, 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 to be resampled.
class RVC_OBJITEM riObjItem;					# objItem for input raster
DlgGetObject("Select raster to resample:", "Raster", riObjItem, "ExistingOnly");


# SET PIPELINE SOURCE FOR INPUT IMAGE
# set input raster as source
class IMAGE_PIPELINE_SOURCE_RVC source_In( riObjItem );
err = source_In.Initialize();
if (err < 0) 
	ReportError(_context.CurrentLineNum, err);
else
	print("Pipeline source intialized.");

printf("Source image has %d lines and %d columns.\n",  source_In.GetTotalRows(), source_In.GetTotalColumns() );


# CHECK THAT SOURCE IMAGE HAS VALID COORDINATE REFERENCE SYSTEM
class IMAGE_PIPELINE_GEOREFERENCE sourceGeoref;
sourceGeoref = source_In.GetGeoreference();

# get coordinate reference system from the source georeference
class SR_COORDREFSYS crsIn;						# input raster's coordinate reference system
crsIn = sourceGeoref.GetCRS();

if (crsIn.IsDefined() == 0 or crsIn.IsLocal() ) {
	PopupMessage("Source coordinate reference system is undefined or local; exiting script.");
	Exit();
	}
else 
	printf("Input image coordinate reference system: %s\n", crsIn.Name );


# CHOOSE REFERENCE RASTER TO RESAMPLE TO
class RVC_OBJITEM refObjItem;
DlgGetObject("Select reference raster:", "Raster", refObjItem, "ExistingOnly");


# SET PIPELINE SOURCE FOR REFERENCE IMAGE
class IMAGE_PIPELINE_SOURCE_RVC source_Ref( refObjItem );
err = source_Ref.Initialize();
if (err < 0) 
	ReportError(_context.CurrentLineNum, err);
else
	print("Reference source initialized.");

printf("Reference image has %d lines and %d columns.\n",  source_Ref.GetTotalRows(), source_Ref.GetTotalColumns() );


# CHECK THAT REFERENCE IMAGE HAS VALID COORDINATE REFERENCE SYSTEM
class IMAGE_PIPELINE_GEOREFERENCE refGeoref;
refGeoref = source_Ref.GetGeoreference();

# get coordinate reference system from the source georeference
class SR_COORDREFSYS crsRef;						# reference raster's coordinate reference system
crsRef = refGeoref.GetCRS();

if (crsRef.IsDefined() == 0 or crsRef.IsLocal() ) {
	PopupMessage("Reference coordinate reference system is undefined or local; exiting script.");
	Exit();
	}
else 
	printf("Reference coordinate reference system: %s\n", crsRef.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 process.
class RVC_OBJITEM rastOutObjItem;
DlgGetObject("Choose raster for resampled output", "Raster", rastOutObjItem, "NewOnly");


# DETERMINE APPROPRIATE RESAMPLING METHOD 
# compare cell sizes (in meters) of input and reference image
class POINT2D scaleIn, scaleRef;		# line and column cell sizes of input and reference as x and y values of POINT2D
class POINT2D locIn, locRef;			# column and line location for which to obtain cell size
locIn.x = source_In.GetTotalColumns() / 2;	# location at center of input image
locIn.y = source_In.GetTotalRows() / 2;
locRef.x = source_Ref.GetTotalColumns() / 2;	# location at center of reference image
locRef.y = source_Ref.GetTotalRows() / 2;

sourceGeoref.ComputeScale(locIn, scaleIn, 1);	# cell sizes of input image in meters
refGeoref.ComputeScale(locRef, scaleRef, 1);		# cell sizes of reference image in meters
scaleIn.y = abs(scaleIn.y);	scaleRef.y = abs(scaleRef.y);

printf("Input image cell sizes: line = %.2f m, col = %.2f m\n", scaleIn.y, scaleIn.x);
printf("Reference image cell sizes: line = %.2f m, col = %.2f m\n", scaleRef.y, scaleRef.x);

numeric inCellArea, outCellArea;					# compute cell areas of input and reference for comparison
inCellArea = scaleIn.x * scaleIn.y;
outCellArea = scaleRef.x * scaleRef.y;

string rsmpMethod$;			# set resampling method
if (outCellArea > 2 * inCellArea) then rsmpMethod$ = "Average";
	else rsmpMethod$ = "CubicConvolution";
printf("Resampling method: %s\n", rsmpMethod$); 


# PIPELINE FILTER to resample source image
class IMAGE_PIPELINE_FILTER_RESAMPLE filter_rsmp(source_In, source_Ref, rsmpMethod$);
err = filter_rsmp.Initialize();
if (err < 0) 
	ReportError(_context.CurrentLineNum, err);
else
	print("Resample filter initialized.");


# PIPELINE TARGET
class IMAGE_PIPELINE_TARGET_RVC target_rvc(filter_rsmp, rastOutObjItem);
target_rvc.SetCompression("DPCM", 0);
err = target_rvc.Initialize();
if (err < 0) 
	ReportError(_context.CurrentLineNum, err);
else
	print("Target initialized.");


# EXECUTE pipeline process
print("Processing...");
target_rvc.Process();

print("Done.");