SetLabelizeParameters1(p1);
// Step 2 default values
- SetUpperThresholdForTrachea(-900);
+ SetTracheaSeedAlgorithm(0);
+ SetUpperThresholdForTrachea(-700);
SetMultiplierForTrachea(5);
SetThresholdStepSizeForTrachea(64);
SetNumberOfSlicesToSkipBeforeSearchingSeed(0);
TracheaVolumeMustBeCheckedFlagOn();
+ SetNumSlices(50);
+ SetMaxElongation(0.5);
+ SetSeedPreProcessingThreshold(-400);
+
// Step 3 default values
SetNumberOfHistogramBins(500);
StartNewStep("Search for the trachea");
SearchForTrachea();
PrintMemory(GetVerboseMemoryFlag(), "After SearchForTrachea");
+ if (m_Seeds.empty()) {
+ clitkExceptionMacro("No seeds for trachea... Aborting.");
+ }
//--------------------------------------------------------------------
//--------------------------------------------------------------------
bool is_orientation_superior(itk::SpatialOrientation::ValidCoordinateOrientationFlags orientation)
{
itk::SpatialOrientation::CoordinateTerms sup = itk::SpatialOrientation::ITK_COORDINATE_Superior;
- std::cout << "orientation: " << std::hex << orientation << "; superior: " << std::hex << sup << std::endl;
- std::cout << std::dec;
-
bool primary = (orientation & 0x0000ff) == sup;
bool secondary = ((orientation & 0x00ff00) >> 8) == sup;
bool tertiary = ((orientation & 0xff0000) >> 16) == sup;
SearchForTracheaSeed2(int numberOfSlices)
{
if (m_Seeds.size() == 0) { // try to find seed (if not zero, it is given by user)
+ if (GetVerboseRegionGrowingFlag())
+ std::cout << "SearchForTracheaSeed2(" << numberOfSlices << ", " << GetMaxElongation() << ")" << std::endl;
+
typedef unsigned char MaskPixelType;
typedef itk::Image<MaskPixelType, ImageType::ImageDimension> MaskImageType;
typedef itk::Image<typename MaskImageType::PixelType, 2> MaskImageType2D;
// threshold to isolate airawys and lungs
typename ThresholdFilterType::Pointer threshold = ThresholdFilterType::New();
threshold->SetLowerThreshold(-2000);
- threshold->SetUpperThreshold(-400);
+ threshold->SetUpperThreshold(GetSeedPreProcessingThreshold());
threshold->SetInput(working_input);
threshold->Update();
writer->SetFileName(file_name.str().c_str());
writer->Update();
}
-
- typename LabelImageType::LabelObjectContainerType shapes_map = label_map->GetLabelObjectContainer();
- typename LabelImageType::LabelObjectContainerType::const_iterator s;
+
typename ShapeLabelType::Pointer shape, max_e_shape;
- double max_elongation = 0.5;
+ double max_elongation = GetMaxElongation();
double max_size = size[0]*size[1]/128;
double max_e = 0;
int nshapes = 0;
- for (s = shapes_map.begin(); s != shapes_map.end(); s++) {
- shape = s->second;
- if (shape->GetSize() > 150 && shape->GetSize() <= max_size) {
+ unsigned int nlables = label_map->GetNumberOfLabelObjects();
+ for (unsigned int j = 0; j < nlables; j++) {
+ shape = label_map->GetNthLabelObject(j);
+ if (shape->Size() > 150 && shape->Size() <= max_size) {
+#if ITK_VERSION_MAJOR < 4
double e = 1 - 1/shape->GetBinaryElongation();
+#else
+ double e = 1 - 1/shape->GetElongation();
+#endif
//double area = 1 - r->Area() ;
if (e < max_elongation) {
nshapes++;
p2[2] = prev_e_centre[2];
double mag = (p2 - p1).GetNorm();
- if (GetVerboseStepFlag()) {
+ if (GetVerboseRegionGrowingFlag()) {
cout.precision(3);
cout << index[2] << ": ";
cout << "region(" << max_e_centre[0] << ", " << max_e_centre[1] << ", " << max_e_centre[2] << "); ";
prev_e_centre= max_e_centre;
}
+ else {
+ if (GetVerboseRegionGrowingFlag()) {
+ cout << "No shapes found at slice " << index[2] << std::endl;
+ }
+ }
}
size_t longest = 0;
}
}
- if (GetVerboseStepFlag())
- std::cout << "seed at: " << trachea_centre << std::endl;
- m_Seeds.push_back(trachea_centre);
+ if (longest > 0) {
+ if (GetVerboseRegionGrowingFlag())
+ std::cout << "seed at: " << trachea_centre << std::endl;
+ m_Seeds.push_back(trachea_centre);
+ }
}
return (m_Seeds.size() != 0);
// compute trachea volume
// if volume not plausible -> skip more slices and restart
+ bool has_seed;
bool stop = false;
double volume = 0.0;
int skip = GetNumberOfSlicesToSkipBeforeSearchingSeed();
while (!stop) {
stop = true;
- if (SearchForTracheaSeed2(50)) {
+ if (GetTracheaSeedAlgorithm() == 0)
+ has_seed = SearchForTracheaSeed(skip);
+ else
+ has_seed = SearchForTracheaSeed2(GetNumSlices());
+
+ if (has_seed) {
TracheaRegionGrowing();
volume = ComputeTracheaVolume()/1000; // assume mm3, so divide by 1000 to get cc
if (GetWriteStepFlag()) {
std::cout << "\t Found trachea with volume " << volume << " cc." << std::endl;
}
}
- else {
- if (GetVerboseStepFlag()) {
- std::cout << "\t The volume of the trachea (" << volume
- << " cc) seems not correct. I skip some slices (" << skip << ") and restart to find seeds."
- << std::endl;
- }
- skip += 5;
- stop = false;
- // empty the list of seed
- m_Seeds.clear();
+ else
+ if (GetTracheaSeedAlgorithm() == 0) {
+ if (GetVerboseStepFlag()) {
+ std::cout << "\t The volume of the trachea (" << volume
+ << " cc) seems not correct. I skip some slices (" << skip << ") and restart to find seeds."
+ << std::endl;
+ }
+ skip += 5;
+ stop = false;
+ // empty the list of seed
+ m_Seeds.clear();
}
if (skip > 0.5 * working_input->GetLargestPossibleRegion().GetSize()[2]) {
// we want to skip more than a half of the image, it is probably a bug