Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ Orientation Analysis (Cleanup)

## Description

This filter will calculate the misorientation between a cell and its 6 face neighbor cells if the *Mask* value is true of the target cell. If the misorientation falls below the user defined threshold for the user defined minimum number of neighbors, then the target cell will have its *Mask* value changed from false to true.
This filter will calculate the misorientation between a cell and each of its 6 face neighbor cells if the *Mask* value of the target cell is **FALSE**. *If the misorientation falls below* the user defined threshold of the *minimum number of neighbors*, then the target cell will have its *Mask* value **CHANGED** from *false* to ***TRUE***.

The filter will iteratively reduce the required number of neighbors from 6 until it reaches the user defined number. So, if the user selects a required number of neighbors of 4, then the filter will run with a required number of neighbors of 6, then 5, then 4 before finishing.

Expand All @@ -16,12 +16,10 @@ Only the *Mask* value defining the cell as *good* or *bad* is changed. No other

## Example Data

| Example Input/Output Images |
|--------------------------------|
| ![](Images/BadDataNeighborOrientationCheckFilter_1.png) |
| The Small IN100 data just after initial alignment filters have completed. |
| ![](Images/BadDataNeighborOrientationCheckFilter_2.png) |
| The Small IN100 data just after running this filter with a *Misorientation Tolerance* of 5 degrees and a *Required Number of Neighbors* of 4. |
| Example Input Image | Example Output Image |
|---------------------|----------------------|
| ![](Images/BadDataNeighborOrientationCheckFilter_1.png) | ![](Images/BadDataNeighborOrientationCheckFilter_2.png) |
| The Small IN100 data just after initial alignment filters have completed. | The Small IN100 data just after running this filter with a *Misorientation Tolerance* of 5 degrees and a *Required Number of Neighbors* of 4. |

From the above before and after images you can see that this filter can help modify a mask that was generated through a simple threshold filter. This filter essentially uses neighbors to determine if a cell point should have had a mask value of false. The majority of cells that were changed from *false* or the black voxels, into valid IPF colored voxels, had a confidence index that fell just below the initial threshold applied (Confidence Index > 0.1 and Image Quality > 120). This filter determines that enough of that cell's neighbors had a mask value of true, the misorientation was < 5 degrees and the cell had enough valid neighbors that the cell's *mask* value was changed from **false** to **true**.

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,12 +2,10 @@

#include "simplnx/Common/Numbers.hpp"
#include "simplnx/DataStructure/DataArray.hpp"
#include "simplnx/DataStructure/DataGroup.hpp"
#include "simplnx/DataStructure/Geometry/ImageGeom.hpp"
#include "simplnx/Utilities/MaskCompareUtilities.hpp"
#include "simplnx/Utilities/MessageHelper.hpp"

#include <EbsdLib/Core/Orientation.hpp>
#include <EbsdLib/LaueOps/LaueOps.h>

using namespace nx::core;
Expand All @@ -34,177 +32,179 @@ const std::atomic_bool& BadDataNeighborOrientationCheck::getCancel()
// -----------------------------------------------------------------------------
Result<> BadDataNeighborOrientationCheck::operator()()
{
float misorientationTolerance = m_InputValues->MisorientationTolerance * numbers::pi_v<float> / 180.0f;
const float misorientationTolerance = m_InputValues->MisorientationTolerance * numbers::pi_v<float> / 180.0f;

auto* imageGeomPtr = m_DataStructure.getDataAs<ImageGeom>(m_InputValues->ImageGeomPath);
const auto* imageGeomPtr = m_DataStructure.getDataAs<ImageGeom>(m_InputValues->ImageGeomPath);
SizeVec3 udims = imageGeomPtr->getDimensions();
const auto& cellPhases = m_DataStructure.getDataRefAs<Int32Array>(m_InputValues->CellPhasesArrayPath);
const auto& quats = m_DataStructure.getDataRefAs<Float32Array>(m_InputValues->QuatsArrayPath);
const auto& crystalStructures = m_DataStructure.getDataRefAs<UInt32Array>(m_InputValues->CrystalStructuresArrayPath);
size_t totalPoints = quats.getNumberOfTuples();
const usize totalPoints = quats.getNumberOfTuples();

std::unique_ptr<MaskCompareUtilities::MaskCompare> maskCompare = nullptr;
std::unique_ptr<MaskCompareUtilities::MaskCompare> maskCompare;
try
{
maskCompare = MaskCompareUtilities::InstantiateMaskCompare(m_DataStructure, m_InputValues->MaskArrayPath);
} catch(const std::out_of_range& exception)
{
// This really should NOT be happening as the path was verified during preflight BUT we may be calling this from
// somewhere else that is NOT going through the normal nx::core::IFilter API of Preflight and Execute
std::string message = fmt::format("Mask Array DataPath does not exist or is not of the correct type (Bool | UInt8) {}", m_InputValues->MaskArrayPath.toString());
return MakeErrorResult(-54900, message);
return MakeErrorResult(-54900, fmt::format("Mask Array DataPath does not exist or is not of the correct type (Bool | UInt8) {}", m_InputValues->MaskArrayPath.toString()));
}

int64_t dims[3] = {
static_cast<int64_t>(udims[0]),
static_cast<int64_t>(udims[1]),
static_cast<int64_t>(udims[2]),
int64 dims[3] = {
static_cast<int64>(udims[0]),
static_cast<int64>(udims[1]),
static_cast<int64>(udims[2]),
};

int32_t good = 1;
int64_t neighbor = 0;
int64_t column = 0, row = 0, plane = 0;
int64 column = 0, row = 0, plane = 0;

int64_t neighpoints[6] = {0, 0, 0, 0, 0, 0};
neighpoints[0] = static_cast<int64_t>(-dims[0] * dims[1]);
neighpoints[1] = static_cast<int64_t>(-dims[0]);
neighpoints[2] = static_cast<int64_t>(-1);
neighpoints[3] = static_cast<int64_t>(1);
neighpoints[4] = static_cast<int64_t>(dims[0]);
neighpoints[5] = static_cast<int64_t>(dims[0] * dims[1]);

float w = 10000.0f;
int64 neighpoints[6] = {0, 0, 0, 0, 0, 0};
neighpoints[0] = -dims[0] * dims[1];
neighpoints[1] = -dims[0];
neighpoints[2] = static_cast<int64>(-1);
neighpoints[3] = static_cast<int64>(1);
neighpoints[4] = dims[0];
neighpoints[5] = dims[0] * dims[1];

std::vector<ebsdlib::LaueOps::Pointer> orientationOps = ebsdlib::LaueOps::GetAllOrientationOps();

std::vector<int32_t> neighborCount(totalPoints, 0);
std::vector<int32> neighborCount(totalPoints, 0);

MessageHelper messageHelper(m_MessageHandler);
ThrottledMessenger throttledMessenger = messageHelper.createThrottledMessenger();
for(size_t i = 0; i < totalPoints; i++)
// Loop over every point finding the number of neighbors that fall within the
// user defined angle tolerance.
for(usize voxelIdx = 0; voxelIdx < totalPoints; voxelIdx++)
{
throttledMessenger.sendThrottledMessage([&]() { return fmt::format("Processing Data {:.2f}% completed", CalculatePercentComplete(i, totalPoints)); });

if(!maskCompare->isTrue(i))
throttledMessenger.sendThrottledMessage([&] { return fmt::format("Processing Data {:.2f}% completed", CalculatePercentComplete(voxelIdx, totalPoints)); });
// If the mask was set to false, then we check this voxel
if(!maskCompare->isTrue(voxelIdx))
{
column = i % dims[0];
row = (i / dims[0]) % dims[1];
plane = i / (dims[0] * dims[1]);
for(int32_t j = 0; j < 6; j++)
// We precalculate the positive voxel quaternion and laue class here to prevent reading and recalculating it for each face below
ebsdlib::QuatD quat1(quats[voxelIdx * 4], quats[voxelIdx * 4 + 1], quats[voxelIdx * 4 + 2], quats[voxelIdx * 4 + 3]);
quat1.positiveOrientation();
const uint32 laueClass1 = crystalStructures[cellPhases[voxelIdx]];

column = voxelIdx % dims[0];
row = (voxelIdx / dims[0]) % dims[1];
plane = voxelIdx / (dims[0] * dims[1]);

// Check the 6 Faces of the voxel
for(int32 faceIdx = 0; faceIdx < 6; faceIdx++)
{
good = 1;
neighbor = i + neighpoints[j];
if(j == 0 && plane == 0)
int64 neighborIdx = static_cast<int64>(voxelIdx) + neighpoints[faceIdx];
// clang-format off
if((faceIdx == 0 && plane == 0) ||
(faceIdx == 1 && row == 0) ||
(faceIdx == 2 && column == 0) ||
(faceIdx == 3 && column == (dims[0] - 1)) ||
(faceIdx == 4 && row == (dims[1] - 1)) ||
(faceIdx == 5 && plane == (dims[2] - 1)))
{
good = 0;
continue;
}
if(j == 5 && plane == (dims[2] - 1))
{
good = 0;
}
if(j == 1 && row == 0)
{
good = 0;
}
if(j == 4 && row == (dims[1] - 1))
{
good = 0;
}
if(j == 2 && column == 0)
{
good = 0;
}
if(j == 3 && column == (dims[0] - 1))
{
good = 0;
}
if(good == 1 && maskCompare->isTrue(neighbor))
{
uint32 laueClass1 = crystalStructures[cellPhases[i]];
ebsdlib::QuatD quat1(quats[i * 4], quats[i * 4 + 1], quats[i * 4 + 2], quats[i * 4 + 3]);
ebsdlib::QuatD quat2(quats[neighbor * 4], quats[neighbor * 4 + 1], quats[neighbor * 4 + 2], quats[neighbor * 4 + 3]);
// clang-format on

if(cellPhases[i] == cellPhases[neighbor] && cellPhases[i] > 0)
// Now compare the mask of the neighbor. If the mask is TRUE, i.e., that voxel
// did not fail the threshold filter that most likely produced the mask array,
// then we can look at that voxel.
if(maskCompare->isTrue(neighborIdx))
{
// Both Cell Phases MUST be the same and be a valid Phase
if(cellPhases[voxelIdx] == cellPhases[neighborIdx] && cellPhases[voxelIdx] > 0)
{
ebsdlib::QuatD quat2(quats[neighborIdx * 4], quats[neighborIdx * 4 + 1], quats[neighborIdx * 4 + 2], quats[neighborIdx * 4 + 3]);
quat2.positiveOrientation();
// Compute the Axis_Angle misorientation between those 2 quaternions
ebsdlib::AxisAngleDType axisAngle = orientationOps[laueClass1]->calculateMisorientation(quat1, quat2);
w = axisAngle[3];
}
if(w < misorientationTolerance)
{
neighborCount[i]++;
// if the angle is less than our tolerance, then we increment the neighbor count
// for this voxel
if(axisAngle[3] < misorientationTolerance)
{
neighborCount[voxelIdx]++;
}
}
}
}
}
}

const int32_t startLevel = 6;
int32_t currentLevel = startLevel;
int32_t counter = 0;
constexpr int32 startLevel = 6;
int32 currentLevel = startLevel;
int32 counter = 0;

while(currentLevel > m_InputValues->NumberOfNeighbors)
// Now we loop over all the points again, but this time we do it as many times
// as the user has requested to iteratively flip voxels
while(currentLevel >= m_InputValues->NumberOfNeighbors)
{
counter = 1;
int32_t loopNumber = 0;
int32 loopNumber = 0;
while(counter > 0)
{
counter = 0;
for(size_t i = 0; i < totalPoints; i++)
counter = 0; // Set this while control variable to zero
for(usize voxelIdx = 0; voxelIdx < totalPoints; voxelIdx++)
{
throttledMessenger.sendThrottledMessage([&]() {
throttledMessenger.sendThrottledMessage([&] {
return fmt::format("Level '{}' of '{}' || Processing Data ('{}') {:.2f}% completed", (startLevel - currentLevel) + 1, startLevel - m_InputValues->NumberOfNeighbors, loopNumber,
CalculatePercentComplete(i, totalPoints));
CalculatePercentComplete(voxelIdx, totalPoints));
});

if(neighborCount[i] >= currentLevel && !maskCompare->isTrue(i))
// We not compare the number-of-neighbors of the current voxel and if it
// is > the current level and the mask is FALSE, then we drop into this
// conditional. The first thing that happens in the conditional is that
// the current voxel's mask value is set to TRUE.
if(neighborCount[voxelIdx] >= currentLevel && !maskCompare->isTrue(voxelIdx))
{
maskCompare->setValue(i, true);
counter++;
column = i % dims[0];
row = (i / dims[0]) % dims[1];
plane = i / (dims[0] * dims[1]);
for(int64_t j = 0; j < 6; j++)
maskCompare->setValue(voxelIdx, true); // current voxel's mask value is set to TRUE.
counter++; // Increment the `counter` to force the loop to iterate again

// We precalculate the positive voxel quaternion and laue class here to prevent reading and recalculating it for each face below
ebsdlib::QuatD quat1(quats[voxelIdx * 4], quats[voxelIdx * 4 + 1], quats[voxelIdx * 4 + 2], quats[voxelIdx * 4 + 3]);
quat1.positiveOrientation();
const uint32 laueClass1 = crystalStructures[cellPhases[voxelIdx]];

// This whole section below is to now look at the neighbor voxels of the
// current voxel that just got flipped to true. This is needed because
// if any of those neighbors mask was `false` then its neighbor count
// is now not correct and will be off-by-one. So we run _almost_ the same
// loop code as above but checking the specific neighbors of the current
// voxel. This part should be termed the "Update Neighbor's Neighbor Count"
column = voxelIdx % dims[0]; // Calculate the column, row, plane
row = (voxelIdx / dims[0]) % dims[1];
plane = voxelIdx / (dims[0] * dims[1]);

for(int64 j = 0; j < 6; j++) // Loop over each of the 6 neighbor faces
{
good = 1;
neighbor = i + neighpoints[j];
if(j == 0 && plane == 0)
{
good = 0;
}
if(j == 5 && plane == (dims[2] - 1))
{
good = 0;
}
if(j == 1 && row == 0)
{
good = 0;
}
if(j == 4 && row == (dims[1] - 1))
int64 neighborIdx = static_cast<int64>(voxelIdx) + neighpoints[j];
// clang-format off
// Do NOT even look at any voxel along the outside boundary of the volume
if((j == 0 && plane == 0) ||
(j == 1 && row == 0) ||
(j == 2 && column == 0) ||
(j == 3 && column == (dims[0] - 1)) ||
(j == 4 && row == (dims[1] - 1)) ||
(j == 5 && plane == (dims[2] - 1)))
{
good = 0;
continue;
}
if(j == 2 && column == 0)
{
good = 0;
}
if(j == 3 && column == (dims[0] - 1))
{
good = 0;
}
if(good == 1 && !maskCompare->isTrue(neighbor))
{
ebsdlib::QuatD quat1(quats[i * 4], quats[i * 4 + 1], quats[i * 4 + 2], quats[i * 4 + 3]);
ebsdlib::QuatD quat2(quats[neighbor * 4], quats[neighbor * 4 + 1], quats[neighbor * 4 + 2], quats[neighbor * 4 + 3]);
// clang-format on

if(cellPhases[i] == cellPhases[neighbor] && cellPhases[i] > 0)
// If the neighbor voxel's mask is false then ....
if(!maskCompare->isTrue(neighborIdx))
{
// Make sure both cell's phase are identical and valid
if(cellPhases[voxelIdx] == cellPhases[neighborIdx] && cellPhases[voxelIdx] > 0)
{
uint32 laueClass1 = crystalStructures[cellPhases[i]];
ebsdlib::QuatD quat2(quats[neighborIdx * 4], quats[neighborIdx * 4 + 1], quats[neighborIdx * 4 + 2], quats[neighborIdx * 4 + 3]);
quat2.positiveOrientation();
// Quaternion Math is not commutative so do not reorder
ebsdlib::AxisAngleDType axisAngle = orientationOps[laueClass1]->calculateMisorientation(quat1, quat2);
w = axisAngle[3];
}
if(w < misorientationTolerance)
{
neighborCount[neighbor]++;
if(axisAngle[3] < misorientationTolerance)
{
neighborCount[neighborIdx]++;
}
}
}
}
Expand Down
Loading
Loading