Thursday, April 10, 2008

Creating a Mosaic of FalconView DTED in the GIS Editor

In my previous post on the subject of DTED in the FalconView ArcGIS Editor, I mentioned that bringing DTED into the GIS Editor is a two step process.  The FalconView user first has the Editor create a feature class that shows all the DTED tiles in the FalconView map data manager; the second step is to select the area of interest and create a DTED mosaic from that area.  The previous post examined the creation of a feature class to show DTED coverage.  This post will examine the technical details of how the GIS Editor uses the geoprocessing capabilities of ArcGIS to add a nicely rendered DTED mosaic to the map.

(Please note that the code below is still somewhat experimental.  Comments and suggestions are welcome.  As FalconView 4.2 is still under development, there are some discussions within Georgia Tech and our user community as to how the final user experience for this will happen.)

1. Get a list of DTED tiles to add to the mosaic

The code block below shows how to enumerate all features the layer of FalconView DTED coverage, collecting a list of DTED tiles which will later be passed to the Mosaic geoprocessing tool.  Note that the only end result from this code block is sMosaicInput, the colon separated list of tiles to add to the mosaic (e.g. "C:\dted\w086\n35.dt1;C:\dted\w086\n34.dt1").  I have highlighted the two most interesting parts of the code: the part that builds the lookup table of selected features and the part that builds sMosaicInput from the features in the lookup table.
string sMosaicInput = null;

// ... some redacted code ...

// The block below populates sMosaicInput

try
{
// Create a lookup table of selected features

IFeatureSelection pFeatureSelection = (IFeatureSelection)pLayer;
Dictionary<int, object> dictSelectedIDs = new Dictionary<int, object>();

IEnumIDs pEnumIDs = pFeatureSelection.SelectionSet.IDs;
for (int id = pEnumIDs.Next(); id != -1; id = pEnumIDs.Next())
dictSelectedIDs[id] = null;

// ... some redacted code ...

// Enumerate all features, adding selected paths to the list

IFeatureLayer pFeatureLayer = (IFeatureLayer)pLayer;
IFeatureCursor pFeatureCursor = pFeatureLayer.FeatureClass.Search(null, true);
int iPathField = pFeatureCursor.FindField("Path");

for (IFeature pFeature = pFeatureCursor.NextFeature(); pFeature != null; pFeature = pFeatureCursor.NextFeature())
{
if (dictSelectedIDs.ContainsKey(pFeature.OID))
{
string sPath = (string)pFeature.get_Value(iPathField);
if (sMosaicInput == null)
sMosaicInput = sPath;
else
sMosaicInput += ';' + sPath;
}
}
}
catch (Exception ex)
{
// ... some redacted code ...

}
2. Use geoprocessing tools to create the mosaic raster

Once we have compiled the string containing the paths of the DTED tiles to be added to the mosaic, creating the mosaic is a simple task requiring the use of two geoprocessing tools.  In the code below, the first highlighted section executes a CreateRasterDataset geoprocessor, creating a scratch raster dataset in a scratch file geodatabase.  The second highlighted lines execute the Mosaic geoprocessor.
// Create a new raster dataset to hold the mosaic

string sRasterName = "dted_mosaic" + Guid.NewGuid().ToString("N");
const string PIXEL_TYPE = "16_BIT_SIGNED"; // Matches DTED 1 & 2 in ArcMap
const int NUMBER_OF_BANDS = 1; // Matches DTED 1 & 2 in ArcMap

CreateRasterDataset createRasterDataset = new CreateRasterDataset(TempGeodatabaseUtils.TempFileGeodatabase, sRasterName, PIXEL_TYPE, NUMBER_OF_BANDS);

Geoprocessor geoprocessor = new Geoprocessor();
geoprocessor.Execute(createRasterDataset, null);
for (int j = 0; j < geoprocessor.MessageCount; j++)
   LoggingUtils.TraceMessage(geoprocessor.GetMessage(j));

// Create the mosaic

Mosaic mosaic = new Mosaic(sMosaicInput, TempGeodatabaseUtils.TempFileGeodatabase + '\\' + sRasterName);
mosaic.mosaic_type = "BLEND";
mosaic.nodata_value = -32767;  // NoData for DTED 1 & 2 in ArcMap

geoprocessor = new Geoprocessor();
geoprocessor.Execute(mosaic, null);
for (int j = 0; j < geoprocessor.MessageCount; j++)
   LoggingUtils.TraceMessage(geoprocessor.GetMessage(j));
3. Symbolize the mosaic and add it to the map

Once the mosaic is created, it must be symbolized before it is added to the FalconView map.  After some experimentation, we found that a simple color ramp shows the features of the DTED well.  (Note that a possible improvement to this would be to add a hillshade effect.)  The first section of highlighted code below shows how to create raster statistics programmatically.  The raster statistics will be used to determine the parameters of the renderer.  The highlighted section shows how to build a simple stretch color ramp.  The last lines show how to create a raster layer from the raw raster, apply the renderer, and add the raster to the map control.
// Set up a raster data object from the mosaic

IWorkspaceFactory pWorkspaceFactory = new FileGDBWorkspaceFactory();
IWorkspace pWorkspace = pWorkspaceFactory.OpenFromFile(TempGeodatabaseUtils.TempFileGeodatabase, hWnd);
IRasterWorkspaceEx pRasterWorkspace = (IRasterWorkspaceEx)pWorkspace;
IRasterDataset pRasterDataset = pRasterWorkspace.OpenRasterDataset(sRasterName);

// Create raster statistics so we can find the highest and lowest points on the raster for the color ramp

IRasterBandCollection pRasterBandCollection = (IRasterBandCollection)pRasterDataset;
IEnumRasterBand pEnumRasterBand = pRasterBandCollection.Bands;
IRasterBand pRasterBand = pEnumRasterBand.Next();
pRasterBand.ComputeStatsAndHist();
IRasterStatistics pRasterStatistics = pRasterBand.Statistics;

// Create a renderer to be used with the raster layer

IRgbColor pFromColor = new RgbColorClass();
pFromColor.Red = 254;
pFromColor.Green = 254;
pFromColor.Blue = 255;

IRgbColor pToColor = new RgbColorClass();
pToColor.Red = 0;
pToColor.Green = 0;
pToColor.Blue = 255;

IAlgorithmicColorRamp pAlgorithmicColorRamp = new AlgorithmicColorRampClass();
pAlgorithmicColorRamp.FromColor = pFromColor;
pAlgorithmicColorRamp.ToColor = pToColor;
pAlgorithmicColorRamp.Algorithm = esriColorRampAlgorithm.esriCIELabAlgorithm;
pAlgorithmicColorRamp.Size = 50;

bool bOK;
pAlgorithmicColorRamp.CreateRamp(out bOK);
LoggingUtils.Assert(bOK, "Failed to create color ramp");

IRasterStretchColorRampRenderer pRasterStretchColorRampRenderer = new RasterStretchColorRampRendererClass();
pRasterStretchColorRampRenderer.BandIndex = 0;
pRasterStretchColorRampRenderer.LabelHigh = pRasterStatistics.Maximum.ToString();
pRasterStretchColorRampRenderer.LabelLow = pRasterStatistics.Minimum.ToString();
pRasterStretchColorRampRenderer.ColorRamp = pAlgorithmicColorRamp;

// Create the raster layer and add it to the map

IRasterLayer pRasterLayer = new RasterLayerClass();
pRasterLayer.CreateFromDataset(pRasterDataset);
pRasterLayer.Name = "FalconView DTED Mosaic";
pRasterLayer.Renderer = (IRasterRenderer)pRasterStretchColorRampRenderer;

pMapControl2.AddLayer(pRasterLayer, 0);







2 comments:

Anonymous said...
This comment has been removed by a blog administrator.
Anonymous said...
This comment has been removed by a blog administrator.