Tutorial 5 - Creating a random data access edge detection algorithm

Files:

• Tutorial/Tutorial5.h
• Tutorial/Tutorial5.cpp

This tutorial will demonstrate how to randomly access data in a raster element and how to create a new raster element while performing an edge detection algorithm. We'll create a plug-in which performs edge detection and returns a RasterElement with the edge detected data.

If you would like to run this tutorial as you follow along with the code, see Running Tutorials With Opticks. If you would like to build and experiment with this tutorial as you follow along with the code, see Building and Running Tutorial plug-ins. If you are running this tutorial in Opticks, you will need to select "Tutorial 5" from the Tutorial toolbar to execute it.

This tutorial is distinct from the others and does not extend the plug-in used in those tutorials.

The Sobel edge detection algorithm can be calculated using convolution. This algorithm will only do edge detection with the first band in the data cube. The algorithm uses two kernels: one to calculate the vertical gradient, and one to calculate the horizontal gradient.

Vertical Gradient:

-1	0	+1
-2	0	+2
-1	0	+1

Horizontal Gradient:

+1	+2	+1
0	0	0
-1	-2	-1

These two results form a vector which describes the edge direction and magnitude at each point in the image. We are not concerned with edge direction so we will display only the absolute magnitude.

Now let's go over the edge detection function.

   template<typename T>
   void edgeDetection(T* pData, DataAccessor pSrcAcc, int row, int col, int rowSize, int colSize)
   {

The first parameter of this function is templated so that the function can use switchOnEncoding() to set its type. This parameter is then used to store the edge detection result.

      int prevCol = std::max(col - 1, 0);
      int prevRow = std::max(row - 1, 0);
      int nextCol = std::min(col + 1, colSize - 1);
      int nextRow = std::min(row + 1, rowSize - 1);

The DataAccessor does not perform any bounds checking. Because of this, we'll use std::max and std::min to cap the convolution filter to the cube boundaries.

      pSrcAcc->toPixel(prevRow, prevCol);
      VERIFYNRV(pSrcAcc.isValid());
      T upperLeftVal = *reinterpret_cast<T*>(pSrcAcc->getColumn());

These lines of code will move the DataAccessor to a pixel located in the upper left corner of the window and get its data. The toPixel() function is used to randomly access data in the cube. The toPixel() function does not perform bound or argument checking. We'll then call getColumn() to get the data from the pixel.

      double gx = -1.0 * upperLeftVal + -2.0 * leftVal + -1.0 * lowerLeftVal + 1.0 * upperRightVal + 2.0 *
         rightVal + 1.0 * lowerRightVal;
      double gy = -1.0 * lowerLeftVal + -2.0 * downVal + -1.0 * lowerRightVal + 1.0 * upperLeftVal + 2.0 *
         upVal + 1.0 * upperRightVal;

Next we'll apply the gradient kernels.

      double magnitude = sqrt(gx * gx + gy * gy);

Now we'll calculate the magnitude of the edge vector.

      *pData = static_cast<T>(magnitude);
   }

Finally, we'll store the magnitude in the result cube.

bool Tutorial5::getOutputSpecification(PlugInArgList*& pOutArgList)
{
   pOutArgList = Service<PlugInManagerServices>()->getPlugInArgList();
   VERIFY(pOutArgList != NULL);
   pOutArgList->addArg<RasterElement>("Result", NULL);
   return true;
}

This tutorial will return a RasterElement output argument which will contain the edge detection data.

   if (pDesc->getDataType() == INT4SCOMPLEX || pDesc->getDataType() == FLT8COMPLEX)
   {
      std::string msg = "Edge detection cannot be performed on complex types.";
      pStep->finalize(Message::Failure, msg);
      if (pProgress != NULL) 
      {
         pProgress->updateProgress(msg, 0, ERRORS);
      }
      return false;

In order to simplify the code for this tutorial, we will not support complex types.

   ModelResource<RasterElement> pResultCube(RasterUtilities::createRasterElement(pCube->getName() +
      "_Edge_Detection_Result", pDesc->getRowCount(), pDesc->getColumnCount(), pDesc->getDataType()));

This line of code will create a new RasterElement which will have the same number of columns and rows as the RasterElement we're edge detecting. The new RasterElement will only need one band since we're only going to edge detect a single band. The ModelResource is an RAII ("resource acquisition is initialization") primitive that manages the life of the new RasterElement.

   FactoryResource<DataRequest> pResultRequest;
   pResultRequest->setWritable(true);
   DataAccessor pDestAcc = pResultCube->getDataAccessor(pResultRequest.release());

We create a DataRequest so that we can get a DataAccessor for the result cube. We call DataRequest::setWritable(true) so that we can perform a write to the DataAccessor we receive, which will allow us to store the edge detection result.

         switchOnEncoding(pDesc->getDataType(), edgeDetection, pDestAcc->getColumn(), pSrcAcc, row, col,
            pDesc->getRowCount(), pDesc->getColumnCount());

These lines of code will take a pixel from the source cube and run the edge detection algorithm to get a magnitude based on the data from the nearby pixels.

         pDestAcc->nextColumn();
      }

      pDestAcc->nextRow();

These lines will move the destination DataAccessor to the next row and column.

   if (!isBatch())
   {
      Service<DesktopServices> pDesktop;

      SpatialDataWindow* pWindow = static_cast<SpatialDataWindow*>(pDesktop->createWindow(pResultCube->getName(),
         SPATIAL_DATA_WINDOW));

      SpatialDataView* pView = (pWindow == NULL) ? NULL : pWindow->getSpatialDataView();
      if (pView == NULL)
      {
         std::string msg = "Unable to create view.";
         pStep->finalize(Message::Failure, msg);
         if (pProgress != NULL) 
         {
            pProgress->updateProgress(msg, 0, ERRORS);
         }
         return false;
      }

      pView->setPrimaryRasterElement(pResultCube.get());
      pView->createLayer(RASTER, pResultCube.get());

Once we get the edge detection results, we'll check to see what mode we are running in. If we're running in interactive mode, we'll create a new SpatialDataWindow that will show the new edge detected cube to the user.

   pOutArgList->setPlugInArgValue("Tutorial5_Result", pResultCube.release());

This line will set the result RasterElement as a output argument. We'll call the release function to transfer ownership from the ModelResource to the plug-in caller.