Revert "ENH: Add SST-zoomed mapping to computeIPFIntensity for properly filled density images"

This reverts commit 238756a.

Author: imikejackson

Source/EbsdLib/LaueOps/LaueOps.cpp

@@ -57,7 +57,6 @@
 #include <canvas_ity.hpp>
 
 #include <algorithm> // for std::max
-#include <array>
 #include <chrono>
 #include <exception>
 #include <iomanip>
@@ -944,21 +943,29 @@ ebsdlib::Rgb LaueOps::generateMisorientationColor(const QuatD& q, const QuatD& r
 }
 
 // -----------------------------------------------------------------------------
-std::array<float, 2> LaueOps::adjustFigureOrigin(std::array<float, 2> figureOrigin, int legendWidth, int legendHeight, const std::vector<float>& margins, float fontPtSize,
-                                                 bool generateEntirePlane) const
+std::array<float, 2> LaueOps::adjustFigureOrigin(
+    std::array<float, 2> figureOrigin,
+    int legendWidth, int legendHeight,
+    const std::vector<float>& margins, float fontPtSize,
+    bool generateEntirePlane) const
 {
   return figureOrigin;
 }
 
 // -----------------------------------------------------------------------------
-UInt8ArrayType::Pointer LaueOps::annotateIPFImage(UInt8ArrayType::Pointer triangleImage, int imageDim, int canvasDim, const std::string& title, bool generateEntirePlane) const
+UInt8ArrayType::Pointer LaueOps::annotateIPFImage(
+    UInt8ArrayType::Pointer triangleImage,
+    int imageDim,
+    int canvasDim,
+    const std::string& title,
+    bool generateEntirePlane) const
 {
   const float fontPtSize = static_cast<float>(canvasDim) / 24.0f;
   const std::vector<float> margins = {
-      fontPtSize * 3,                       // Top
-      static_cast<float>(canvasDim / 7.0f), // Right
-      fontPtSize * 2,                       // Bottom
-      static_cast<float>(canvasDim / 7.0f)  // Left
+      fontPtSize * 3,                        // Top
+      static_cast<float>(canvasDim / 7.0f),  // Right
+      fontPtSize * 2,                        // Bottom
+      static_cast<float>(canvasDim / 7.0f)   // Left
   };
 
   int legendHeight = canvasDim - static_cast<int>(margins[0]) - static_cast<int>(margins[2]);
@@ -1006,7 +1013,10 @@ UInt8ArrayType::Pointer LaueOps::annotateIPFImage(UInt8ArrayType::Pointer triang
   context.fill();
 
   // Draw the triangle image onto the canvas
-  context.draw_image(image->getPointer(0), imageDim, imageDim, imageDim * image->getNumberOfComponents(), figureOrigin[0], figureOrigin[1], static_cast<float>(legendWidth),
+  context.draw_image(image->getPointer(0), imageDim, imageDim,
+                     imageDim * image->getNumberOfComponents(),
+                     figureOrigin[0], figureOrigin[1],
+                     static_cast<float>(legendWidth),
                      static_cast<float>(legendHeight));
 
   // Draw title
@@ -1018,14 +1028,20 @@ UInt8ArrayType::Pointer LaueOps::annotateIPFImage(UInt8ArrayType::Pointer triang
   drawIPFAnnotations(context, canvasDim, fontPtSize, margins, figureOrigin, figureCenter, generateEntirePlane);
 
   // Extract rendered pixels and remove alpha channel
-  ebsdlib::UInt8ArrayType::Pointer rgbaCanvasImage = ebsdlib::UInt8ArrayType::CreateArray(canvasDim * canvasDim, {4ULL}, "Annotated IPF", true);
+  ebsdlib::UInt8ArrayType::Pointer rgbaCanvasImage = ebsdlib::UInt8ArrayType::CreateArray(
+      canvasDim * canvasDim, {4ULL}, "Annotated IPF", true);
   context.get_image_data(rgbaCanvasImage->getPointer(0), canvasDim, canvasDim, canvasDim * 4, 0, 0);
 
   return ebsdlib::RemoveAlphaChannel(rgbaCanvasImage.get());
 }
 
 // -----------------------------------------------------------------------------
-UInt8ArrayType::Pointer LaueOps::drawColorBar(UInt8ArrayType::Pointer image, int canvasDim, int numColors, double minValue, double maxValue, bool isMRD) const
+UInt8ArrayType::Pointer LaueOps::drawColorBar(
+    UInt8ArrayType::Pointer image,
+    int canvasDim,
+    int numColors,
+    double minValue, double maxValue,
+    bool isMRD) const
 {
   const float fontPtSize = static_cast<float>(canvasDim) / 24.0f;
 
@@ -1048,7 +1064,10 @@ UInt8ArrayType::Pointer LaueOps::drawColorBar(UInt8ArrayType::Pointer image, int
 
   canvas_ity::canvas context(canvasDim, canvasDim);
   // Put the existing image onto the canvas
-  context.draw_image(rgbaImage->getPointer(0), canvasDim, canvasDim, canvasDim * 4, 0.0f, 0.0f, static_cast<float>(canvasDim), static_cast<float>(canvasDim));
+  context.draw_image(rgbaImage->getPointer(0), canvasDim, canvasDim,
+                     canvasDim * 4, 0.0f, 0.0f,
+                     static_cast<float>(canvasDim),
+                     static_cast<float>(canvasDim));
 
   // Color bar dimensions
   const float barLeft = static_cast<float>(canvasDim) * 0.80f;
@@ -1118,7 +1137,9 @@ UInt8ArrayType::Pointer LaueOps::drawColorBar(UInt8ArrayType::Pointer image, int
 }
 
 // -----------------------------------------------------------------------------
-std::vector<UInt8ArrayType::Pointer> LaueOps::generateAnnotatedIPFDensity(InversePoleFigureConfiguration_t& config, std::pair<double, double>* outMinMax) const
+std::vector<UInt8ArrayType::Pointer> LaueOps::generateAnnotatedIPFDensity(
+    InversePoleFigureConfiguration_t& config,
+    std::pair<double, double>* outMinMax) const
 {
   // Validate square images
   if(config.imageWidth != config.imageHeight)
@@ -1151,29 +1172,10 @@ std::vector<UInt8ArrayType::Pointer> LaueOps::generateAnnotatedIPFDensity(Invers
   ebsdlib::FloatArrayType::Pointer dirs1 = InversePoleFigureUtilities::computeIPFDirections(*this, config.eulers, config.sampleDirections[1]);
   ebsdlib::FloatArrayType::Pointer dirs2 = InversePoleFigureUtilities::computeIPFDirections(*this, config.eulers, config.sampleDirections[2]);
 
-  // Step 2: Compute SST bounding box for zoomed density images
-  // getIpfColorAngleLimits returns {etaMin, etaMax, chiMax(eta)} in radians.
-  auto angleLimits = getIpfColorAngleLimits(0.0);
-  double etaMin = angleLimits[0];
-  double etaMax = angleLimits[1];
-
-  // Get chiMax at etaMax (gives the largest chiMax for cubic; same for others)
-  auto angleLimitsAtMax = getIpfColorAngleLimits(etaMax);
-  double chiMax = angleLimitsAtMax[2];
-
-  // Also check chiMax at etaMin in case it's larger
-  auto angleLimitsAtMin = getIpfColorAngleLimits(etaMin);
-  if(angleLimitsAtMin[2] > chiMax)
-  {
-    chiMax = angleLimitsAtMin[2];
-  }
-
-  std::array<double, 4> sstBBox = {etaMin, etaMax, 0.0, chiMax};
-
-  // Compute intensity images with SST-zoomed mapping
-  ebsdlib::DoubleArrayType::Pointer intensity0 = InversePoleFigureUtilities::computeIPFIntensity(*this, dirs0.get(), imageDim, imageDim, config.lambertDim, config.normalizeMRD, &sstBBox);
-  ebsdlib::DoubleArrayType::Pointer intensity1 = InversePoleFigureUtilities::computeIPFIntensity(*this, dirs1.get(), imageDim, imageDim, config.lambertDim, config.normalizeMRD, &sstBBox);
-  ebsdlib::DoubleArrayType::Pointer intensity2 = InversePoleFigureUtilities::computeIPFIntensity(*this, dirs2.get(), imageDim, imageDim, config.lambertDim, config.normalizeMRD, &sstBBox);
+  // Step 2: Compute intensity images
+  ebsdlib::DoubleArrayType::Pointer intensity0 = InversePoleFigureUtilities::computeIPFIntensity(*this, dirs0.get(), imageDim, imageDim, config.lambertDim, config.normalizeMRD);
+  ebsdlib::DoubleArrayType::Pointer intensity1 = InversePoleFigureUtilities::computeIPFIntensity(*this, dirs1.get(), imageDim, imageDim, config.lambertDim, config.normalizeMRD);
+  ebsdlib::DoubleArrayType::Pointer intensity2 = InversePoleFigureUtilities::computeIPFIntensity(*this, dirs2.get(), imageDim, imageDim, config.lambertDim, config.normalizeMRD);
 
   // Step 3: Find global min/max
   double globalMax = std::numeric_limits<double>::lowest();

Source/EbsdLib/Utilities/InversePoleFigureUtilities.cpp

@@ -137,7 +137,7 @@ ebsdlib::FloatArrayType::Pointer InversePoleFigureUtilities::computeIPFDirection
 
 // -----------------------------------------------------------------------------
 ebsdlib::DoubleArrayType::Pointer InversePoleFigureUtilities::computeIPFIntensity(const LaueOps& ops, ebsdlib::FloatArrayType* ipfDirections, int imageWidth, int imageHeight, int lambertDim,
-                                                                                  bool normalizeMRD, const std::array<double, 4>* sstBoundingBox)
+                                                                                  bool normalizeMRD)
 {
   // Step 1: Bin the crystal directions into the Lambert projection
   float sphereRadius = 1.0f;
@@ -169,112 +169,66 @@ ebsdlib::DoubleArrayType::Pointer InversePoleFigureUtilities::computeIPFIntensit
     // If not MRD, leave as raw counts
   }
 
-  // Step 3: Create the output intensity image
+  // Step 3: Create the output intensity image using equal-area projection
   std::vector<size_t> tDims = {static_cast<size_t>(imageWidth * imageHeight)};
   std::vector<size_t> cDims = {1};
   ebsdlib::DoubleArrayType::Pointer intensity = ebsdlib::DoubleArrayType::CreateArray(tDims, cDims, "IPF_Intensity", true);
   double* intensityPtr = intensity->getPointer(0);
 
-  if(sstBoundingBox != nullptr)
-  {
-    // SST-zoomed mode: map pixels directly to (eta, chi) within the bounding box
-    double etaMin = (*sstBoundingBox)[0];
-    double etaMax = (*sstBoundingBox)[1];
-    double chiMin = (*sstBoundingBox)[2];
-    double chiMax = (*sstBoundingBox)[3];
+  // Lambert azimuthal equal-area projection centered on north pole
+  // Maps the upper hemisphere (z >= 0) to a disk of radius sqrt(2)
+  float unitRadius = std::sqrt(2.0f);
+  float span = 2.0f * unitRadius;
+  float xres = span / static_cast<float>(imageWidth);
+  float yres = span / static_cast<float>(imageHeight);
+
+  int halfWidth = imageWidth / 2;
+  int halfHeight = imageHeight / 2;
 
-    for(int y = 0; y < imageHeight; y++)
+  for(int y = 0; y < imageHeight; y++)
+  {
+    for(int x = 0; x < imageWidth; x++)
     {
-      for(int x = 0; x < imageWidth; x++)
+      int index = y * imageWidth + x;
+
+      // Map pixel to equal-area projection coordinates
+      float xtmp = static_cast<float>(x - halfWidth) * xres + (xres * 0.5f);
+      float ytmp = static_cast<float>(y - halfHeight) * yres + (yres * 0.5f);
+
+      float rhoSq = xtmp * xtmp + ytmp * ytmp;
+
+      // Check if within hemisphere disk
+      if(rhoSq > 2.0f)
       {
-        int index = y * imageWidth + x;
-
-        // Map pixel to (eta, chi) within the bounding box
-        // x maps to eta (left=etaMin, right=etaMax)
-        // y maps to chi (top=chiMin, bottom=chiMax)
-        double eta = etaMin + (static_cast<double>(x) + 0.5) / static_cast<double>(imageWidth) * (etaMax - etaMin);
-        double chi = chiMin + (static_cast<double>(y) + 0.5) / static_cast<double>(imageHeight) * (chiMax - chiMin);
-
-        if(!ops.inUnitTriangle(eta, chi))
-        {
-          intensityPtr[index] = -1.0;
-          continue;
-        }
-
-        // Convert (eta, chi) to unit sphere xyz
-        double sinChi = std::sin(chi);
-        std::array<float, 3> xyz = {static_cast<float>(sinChi * std::cos(eta)), static_cast<float>(sinChi * std::sin(eta)), static_cast<float>(std::cos(chi))};
-
-        // Look up intensity from Lambert bins
-        std::array<float, 2> sqCoord = {0.0f, 0.0f};
-        bool isNorth = lambert->getSquareCoord(xyz.data(), sqCoord.data());
-        if(isNorth)
-        {
-          intensityPtr[index] = lambert->getInterpolatedValue(ModifiedLambertProjection::NorthSquare, sqCoord.data());
-        }
-        else
-        {
-          intensityPtr[index] = lambert->getInterpolatedValue(ModifiedLambertProjection::SouthSquare, sqCoord.data());
-        }
+        intensityPtr[index] = -1.0; // Outside hemisphere
+        continue;
       }
-    }
-  }
-  else
-  {
-    // Default full Lambert hemisphere disk mapping
-    float unitRadius = std::sqrt(2.0f);
-    float span = 2.0f * unitRadius;
-    float xres = span / static_cast<float>(imageWidth);
-    float yres = span / static_cast<float>(imageHeight);
 
-    int halfWidth = imageWidth / 2;
-    int halfHeight = imageHeight / 2;
+      // Inverse Lambert azimuthal equal-area projection (north pole centered)
+      float t = std::sqrt(1.0f - rhoSq / 4.0f);
+      std::array<float, 3> xyz = {xtmp * t, ytmp * t, 1.0f - rhoSq / 2.0f};
 
-    for(int y = 0; y < imageHeight; y++)
-    {
-      for(int x = 0; x < imageWidth; x++)
+      // Compute chi (polar angle from z-axis) and eta (azimuthal angle)
+      double chi = std::acos(static_cast<double>(xyz[2]));
+      double eta = std::atan2(static_cast<double>(xyz[1]), static_cast<double>(xyz[0]));
+
+      // Check if direction is inside the Standard Stereographic Triangle
+      if(!ops.inUnitTriangle(eta, chi))
+      {
+        intensityPtr[index] = -1.0; // Outside SST
+        continue;
+      }
+
+      // Look up the interpolated intensity from the Lambert projection
+      std::array<float, 2> sqCoord = {0.0f, 0.0f};
+      bool isNorth = lambert->getSquareCoord(xyz.data(), sqCoord.data());
+      if(isNorth)
+      {
+        intensityPtr[index] = lambert->getInterpolatedValue(ModifiedLambertProjection::NorthSquare, sqCoord.data());
+      }
+      else
       {
-        int index = y * imageWidth + x;
-
-        // Map pixel to equal-area projection coordinates
-        float xtmp = static_cast<float>(x - halfWidth) * xres + (xres * 0.5f);
-        float ytmp = static_cast<float>(y - halfHeight) * yres + (yres * 0.5f);
-
-        float rhoSq = xtmp * xtmp + ytmp * ytmp;
-
-        // Check if within hemisphere disk
-        if(rhoSq > 2.0f)
-        {
-          intensityPtr[index] = -1.0; // Outside hemisphere
-          continue;
-        }
-
-        // Inverse Lambert azimuthal equal-area projection (north pole centered)
-        float t = std::sqrt(1.0f - rhoSq / 4.0f);
-        std::array<float, 3> xyz = {xtmp * t, ytmp * t, 1.0f - rhoSq / 2.0f};
-
-        // Compute chi (polar angle from z-axis) and eta (azimuthal angle)
-        double chi = std::acos(static_cast<double>(xyz[2]));
-        double eta = std::atan2(static_cast<double>(xyz[1]), static_cast<double>(xyz[0]));
-
-        // Check if direction is inside the Standard Stereographic Triangle
-        if(!ops.inUnitTriangle(eta, chi))
-        {
-          intensityPtr[index] = -1.0; // Outside SST
-          continue;
-        }
-
-        // Look up the interpolated intensity from the Lambert projection
-        std::array<float, 2> sqCoord = {0.0f, 0.0f};
-        bool isNorth = lambert->getSquareCoord(xyz.data(), sqCoord.data());
-        if(isNorth)
-        {
-          intensityPtr[index] = lambert->getInterpolatedValue(ModifiedLambertProjection::NorthSquare, sqCoord.data());
-        }
-        else
-        {
-          intensityPtr[index] = lambert->getInterpolatedValue(ModifiedLambertProjection::SouthSquare, sqCoord.data());
-        }
+        intensityPtr[index] = lambert->getInterpolatedValue(ModifiedLambertProjection::SouthSquare, sqCoord.data());
       }
     }
   }
@@ -283,8 +237,7 @@ ebsdlib::DoubleArrayType::Pointer InversePoleFigureUtilities::computeIPFIntensit
 }
 
 // -----------------------------------------------------------------------------
-void InversePoleFigureUtilities::createIPFColorImage(ebsdlib::DoubleArrayType* intensity, int imageWidth, int imageHeight, int numColors, double minScale, double maxScale,
-                                                     ebsdlib::UInt8ArrayType* rgba)
+void InversePoleFigureUtilities::createIPFColorImage(ebsdlib::DoubleArrayType* intensity, int imageWidth, int imageHeight, int numColors, double minScale, double maxScale, ebsdlib::UInt8ArrayType* rgba)
 {
   // Initialize the image with all zeros
   rgba->initializeWithZeros();

Source/EbsdLib/Utilities/InversePoleFigureUtilities.h

@@ -56,17 +56,17 @@ class LaueOps; // Forward declaration
  */
 struct InversePoleFigureConfiguration_t
 {
-  ebsdlib::FloatArrayType* eulers;            ///<* The Euler Angles (in Radians) to use for the inverse pole figure
-  std::array<Matrix3X1D, 3> sampleDirections; ///<* 3 orthogonal sample reference directions (e.g., RD, TD, ND)
-  int imageWidth;                             ///<* The width of the generated inverse pole figure image in pixels
-  int imageHeight;                            ///<* The height of the generated inverse pole figure image in pixels
-  int lambertDim;                             ///<* The dimensions in voxels of the Lambert Square used for binning/smoothing
-  int numColors;                              ///<* The number of colors to use in the color map
-  std::string colorMap;                       ///<* Name of the ColorMap to use
-  bool normalizeMRD;                          ///<* true=normalize to MRD (Multiples of Random Distribution), false=raw counts
-  std::vector<std::string> labels;            ///<* The labels for each of the 3 inverse pole figures (e.g., "RD", "TD", "ND")
-  std::string phaseName;                      ///<* The name of the phase
-  bool FlipFinalImage;                        ///<* If TRUE, the final image will be flipped across the X Axis so that +Y axis points UP
+  ebsdlib::FloatArrayType* eulers;              ///<* The Euler Angles (in Radians) to use for the inverse pole figure
+  std::array<Matrix3X1D, 3> sampleDirections;   ///<* 3 orthogonal sample reference directions (e.g., RD, TD, ND)
+  int imageWidth;                                ///<* The width of the generated inverse pole figure image in pixels
+  int imageHeight;                               ///<* The height of the generated inverse pole figure image in pixels
+  int lambertDim;                                ///<* The dimensions in voxels of the Lambert Square used for binning/smoothing
+  int numColors;                                 ///<* The number of colors to use in the color map
+  std::string colorMap;                          ///<* Name of the ColorMap to use
+  bool normalizeMRD;                             ///<* true=normalize to MRD (Multiples of Random Distribution), false=raw counts
+  std::vector<std::string> labels;               ///<* The labels for each of the 3 inverse pole figures (e.g., "RD", "TD", "ND")
+  std::string phaseName;                         ///<* The name of the phase
+  bool FlipFinalImage;                           ///<* If TRUE, the final image will be flipped across the X Axis so that +Y axis points UP
 };
 
 /**
@@ -104,13 +104,9 @@ class EbsdLib_EXPORT InversePoleFigureUtilities
    * @param imageHeight Output image height in pixels
    * @param lambertDim Lambert square dimension for binning/smoothing
    * @param normalizeMRD true to normalize to MRD, false for raw counts
-   * @param sstBoundingBox Optional SST bounding box {etaMin, etaMax, chiMin, chiMax} in radians.
-   *   When provided, pixels map to only this region in (eta, chi) space, making the SST fill the image.
-   *   When nullptr, uses the default full Lambert hemisphere disk mapping.
    * @return DoubleArrayType intensity image (imageWidth * imageHeight). Pixels outside SST have value -1.0.
    */
-  static ebsdlib::DoubleArrayType::Pointer computeIPFIntensity(const LaueOps& ops, ebsdlib::FloatArrayType* ipfDirections, int imageWidth, int imageHeight, int lambertDim, bool normalizeMRD,
-                                                                const std::array<double, 4>* sstBoundingBox = nullptr);
+  static ebsdlib::DoubleArrayType::Pointer computeIPFIntensity(const LaueOps& ops, ebsdlib::FloatArrayType* ipfDirections, int imageWidth, int imageHeight, int lambertDim, bool normalizeMRD);
 
   /**
    * @brief Converts an intensity image to RGBA with SST masking. Pixels inside the SST