diff --git a/Modules/DiffusionCore/IODataStructures/mitkFiberBundle.cpp b/Modules/DiffusionCore/IODataStructures/mitkFiberBundle.cpp index 30b7de5..70f3858 100644 --- a/Modules/DiffusionCore/IODataStructures/mitkFiberBundle.cpp +++ b/Modules/DiffusionCore/IODataStructures/mitkFiberBundle.cpp @@ -1,2903 +1,2902 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #include "mitkFiberBundle.h" #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include const char* mitk::FiberBundle::FIBER_ID_ARRAY = "Fiber_IDs"; mitk::FiberBundle::FiberBundle( vtkPolyData* fiberPolyData ) : m_NumFibers(0) { m_TrackVisHeader.hdr_size = 0; m_FiberWeights = vtkSmartPointer::New(); m_FiberWeights->SetName("FIBER_WEIGHTS"); m_FiberPolyData = vtkSmartPointer::New(); if (fiberPolyData != nullptr) m_FiberPolyData = fiberPolyData; else { this->m_FiberPolyData->SetPoints(vtkSmartPointer::New()); this->m_FiberPolyData->SetLines(vtkSmartPointer::New()); } this->UpdateFiberGeometry(); this->GenerateFiberIds(); this->ColorFibersByOrientation(); } mitk::FiberBundle::~FiberBundle() { } mitk::FiberBundle::Pointer mitk::FiberBundle::GetDeepCopy() { mitk::FiberBundle::Pointer newFib = mitk::FiberBundle::New(m_FiberPolyData); newFib->SetFiberColors(this->m_FiberColors); newFib->SetFiberWeights(this->m_FiberWeights); newFib->SetTrackVisHeader(this->GetTrackVisHeader()); return newFib; } vtkSmartPointer mitk::FiberBundle::GeneratePolyDataByIds(std::vector fiberIds, vtkSmartPointer weights) { vtkSmartPointer newFiberPolyData = vtkSmartPointer::New(); vtkSmartPointer newLineSet = vtkSmartPointer::New(); vtkSmartPointer newPointSet = vtkSmartPointer::New(); weights->SetNumberOfValues(fiberIds.size()); int counter = 0; auto finIt = fiberIds.begin(); while ( finIt != fiberIds.end() ) { if (*finIt>GetNumFibers()){ MITK_INFO << "FiberID can not be negative or >NumFibers!!! check id Extraction!" << *finIt; break; } vtkSmartPointer fiber = m_FiberIdDataSet->GetCell(*finIt);//->DeepCopy(fiber); vtkSmartPointer fibPoints = fiber->GetPoints(); vtkSmartPointer newFiber = vtkSmartPointer::New(); newFiber->GetPointIds()->SetNumberOfIds( fibPoints->GetNumberOfPoints() ); for(int i=0; iGetNumberOfPoints(); i++) { newFiber->GetPointIds()->SetId(i, newPointSet->GetNumberOfPoints()); newPointSet->InsertNextPoint(fibPoints->GetPoint(i)[0], fibPoints->GetPoint(i)[1], fibPoints->GetPoint(i)[2]); } weights->InsertValue(counter, this->GetFiberWeight(*finIt)); newLineSet->InsertNextCell(newFiber); ++finIt; ++counter; } newFiberPolyData->SetPoints(newPointSet); newFiberPolyData->SetLines(newLineSet); return newFiberPolyData; } // merge two fiber bundles mitk::FiberBundle::Pointer mitk::FiberBundle::AddBundles(std::vector< mitk::FiberBundle::Pointer > fibs) { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); // add current fiber bundle vtkSmartPointer weights = vtkSmartPointer::New(); auto num_weights = this->GetNumFibers(); for (auto fib : fibs) num_weights += fib->GetNumFibers(); weights->SetNumberOfValues(num_weights); unsigned int counter = 0; for (unsigned int i=0; iGetNumberOfCells(); ++i) { vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights->InsertValue(counter, this->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } for (auto fib : fibs) { // add new fiber bundle for (unsigned int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = fib->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } } // initialize PolyData vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // initialize fiber bundle mitk::FiberBundle::Pointer newFib = mitk::FiberBundle::New(vNewPolyData); newFib->SetFiberWeights(weights); return newFib; } // merge two fiber bundles mitk::FiberBundle::Pointer mitk::FiberBundle::AddBundle(mitk::FiberBundle* fib) { if (fib==nullptr) return this->GetDeepCopy(); MITK_INFO << "Adding fibers"; vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); // add current fiber bundle vtkSmartPointer weights = vtkSmartPointer::New(); weights->SetNumberOfValues(this->GetNumFibers()+fib->GetNumFibers()); unsigned int counter = 0; for (unsigned int i=0; iGetNumberOfCells(); i++) { vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights->InsertValue(counter, this->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } // add new fiber bundle for (unsigned int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = fib->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } // initialize PolyData vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // initialize fiber bundle mitk::FiberBundle::Pointer newFib = mitk::FiberBundle::New(vNewPolyData); newFib->SetFiberWeights(weights); return newFib; } // Only retain fibers with a weight larger than the specified threshold mitk::FiberBundle::Pointer mitk::FiberBundle::FilterByWeights(float weight_thr, bool invert) { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); std::vector weights; for (unsigned int i=0; iGetNumFibers(); i++) { if ( (invert && this->GetFiberWeight(i)>weight_thr) || (!invert && this->GetFiberWeight(i)<=weight_thr)) continue; vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); weights.push_back(this->GetFiberWeight(i)); } // initialize PolyData vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // initialize fiber bundle mitk::FiberBundle::Pointer newFib = mitk::FiberBundle::New(vNewPolyData); for (unsigned int i=0; iSetFiberWeight(i, weights.at(i)); newFib->SetTrackVisHeader(this->GetTrackVisHeader()); return newFib; } // Only retain a subsample of the fibers mitk::FiberBundle::Pointer mitk::FiberBundle::SubsampleFibers(float factor, bool random_seed) { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int new_num_fibs = static_cast(std::round(this->GetNumFibers()*factor)); MITK_INFO << "Subsampling fibers with factor " << factor << "(" << new_num_fibs << "/" << this->GetNumFibers() << ")"; // add current fiber bundle vtkSmartPointer weights = vtkSmartPointer::New(); weights->SetNumberOfValues(new_num_fibs); std::vector< unsigned int > ids; for (unsigned int i=0; iGetNumFibers(); i++) ids.push_back(i); if (random_seed) std::srand(static_cast(std::time(nullptr))); else std::srand(0); std::random_shuffle(ids.begin(), ids.end()); unsigned int counter = 0; for (unsigned int i=0; iGetCell(ids.at(i)); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights->InsertValue(counter, this->GetFiberWeight(ids.at(i))); vNewLines->InsertNextCell(container); counter++; } // initialize PolyData vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // initialize fiber bundle mitk::FiberBundle::Pointer newFib = mitk::FiberBundle::New(vNewPolyData); newFib->SetFiberWeights(weights); newFib->SetTrackVisHeader(this->GetTrackVisHeader()); return newFib; } // subtract two fiber bundles mitk::FiberBundle::Pointer mitk::FiberBundle::SubtractBundle(mitk::FiberBundle* fib) { if (fib==nullptr) return this->GetDeepCopy(); MITK_INFO << "Subtracting fibers"; vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); std::vector< std::vector< itk::Point > > points1; for(unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); if (points==nullptr || numPoints<=0) continue; itk::Point start = mitk::imv::GetItkPoint(points->GetPoint(0)); itk::Point end = mitk::imv::GetItkPoint(points->GetPoint(numPoints-1)); points1.push_back( {start, end} ); } std::vector< std::vector< itk::Point > > points2; for(unsigned int i=0; iGetNumFibers(); i++ ) { vtkCell* cell = fib->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); if (points==nullptr || numPoints<=0) continue; itk::Point start =mitk::imv::GetItkPoint(points->GetPoint(0)); itk::Point end =mitk::imv::GetItkPoint(points->GetPoint(numPoints-1)); points2.push_back( {start, end} ); } // int progress = 0; std::vector< int > ids; #pragma omp parallel for for (int i=0; i(points1.size()); i++) { bool match = false; for (unsigned int j=0; j(i)); auto v2 = points2.at(j); float dist=0; for (unsigned int c=0; cGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); if (points==nullptr || numPoints<=0) continue; vtkSmartPointer container = vtkSmartPointer::New(); for( int j=0; jInsertNextPoint(points->GetPoint(j)); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); } if(vNewLines->GetNumberOfCells()==0) return mitk::FiberBundle::New(); // initialize PolyData vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // initialize fiber bundle return mitk::FiberBundle::New(vNewPolyData); } /* * set PolyData (additional flag to recompute fiber geometry, default = true) */ void mitk::FiberBundle::SetFiberPolyData(vtkSmartPointer fiberPD, bool updateGeometry) { if (fiberPD == nullptr) this->m_FiberPolyData = vtkSmartPointer::New(); else { m_FiberPolyData->CopyStructure(fiberPD); // m_FiberPolyData->DeepCopy(fiberPD); } m_NumFibers = static_cast(m_FiberPolyData->GetNumberOfLines()); if (updateGeometry) UpdateFiberGeometry(); GenerateFiberIds(); ColorFibersByOrientation(); } /* * return vtkPolyData */ vtkSmartPointer mitk::FiberBundle::GetFiberPolyData() const { return m_FiberPolyData; } void mitk::FiberBundle::ColorFibersByLength(bool opacity, bool weight_fibers, mitk::LookupTable::LookupTableType type) { if (m_MaxFiberLength<=0) return; auto numOfPoints = this->GetNumberOfPoints(); //colors and alpha value for each single point, RGBA = 4 components unsigned char rgba[4] = {0,0,0,0}; m_FiberColors = vtkSmartPointer::New(); m_FiberColors->Allocate(numOfPoints * 4); m_FiberColors->SetNumberOfComponents(4); m_FiberColors->SetName("FIBER_COLORS"); auto numOfFibers = m_FiberPolyData->GetNumberOfLines(); if (numOfFibers < 1) return; mitk::LookupTable::Pointer mitkLookup = mitk::LookupTable::New(); mitkLookup->SetType(type); if (type!=mitk::LookupTable::MULTILABEL) mitkLookup->GetVtkLookupTable()->SetTableRange(m_MinFiberLength, m_MaxFiberLength); unsigned int count = 0; for (unsigned int i=0; iGetNumberOfCells(); i++) { vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); float l = m_FiberLengths.at(i)/m_MaxFiberLength; double color[3]; mitkLookup->GetColor(m_FiberLengths.at(i), color); for (int j=0; j(255.0 * color[0]); rgba[1] = static_cast(255.0 * color[1]); rgba[2] = static_cast(255.0 * color[2]); if (opacity) rgba[3] = static_cast(255.0f * l); else rgba[3] = static_cast(255.0); m_FiberColors->InsertTypedTuple(cell->GetPointId(j), rgba); count++; } if (weight_fibers) this->SetFiberWeight(i, m_FiberLengths.at(i)); } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::ColorSinglePoint(int f_idx, int p_idx, double rgb[3]) { // vtkPoints* extrPoints = m_FiberPolyData->GetPoints(); // vtkIdType numOfPoints = 0; // if (extrPoints!=nullptr) // numOfPoints = extrPoints->GetNumberOfPoints(); // //colors and alpha value for each single point, RGBA = 4 components unsigned char rgba[4] = {0,0,0,0}; // m_FiberColors = vtkSmartPointer::New(); // m_FiberColors->Allocate(numOfPoints * 4); // m_FiberColors->SetNumberOfComponents(4); // m_FiberColors->SetName("FIBER_COLORS"); // auto numOfFibers = m_FiberPolyData->GetNumberOfLines(); // if (numOfFibers < 1) // return; // /* extract single fibers of fiberBundle */ // vtkCellArray* fiberList = m_FiberPolyData->GetLines(); // fiberList->InitTraversal(); // for (int fi=0; fiGetNextCell(num_points, idList); // fiberList->GetCell(f_idx, num_points, idList); vtkCell* cell = m_FiberPolyData->GetCell(f_idx); // /* single fiber checkpoints: is number of points valid */ // if (p_idx < num_points) // { rgba[0] = static_cast(255.0 * rgb[0]); rgba[1] = static_cast(255.0 * rgb[1]); rgba[2] = static_cast(255.0 * rgb[2]); rgba[3] = 255; m_FiberColors->InsertTypedTuple(cell->GetPointId(p_idx), rgba); // } // } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::ColorFibersByOrientation() { //===== FOR WRITING A TEST ======================== // colorT size == tupelComponents * tupelElements // compare color results // to cover this code 100% also PolyData needed, where colorarray already exists // + one fiber with exactly 1 point // + one fiber with 0 points //================================================= vtkPoints* extrPoints = m_FiberPolyData->GetPoints(); vtkIdType numOfPoints = 0; if (extrPoints!=nullptr) numOfPoints = extrPoints->GetNumberOfPoints(); //colors and alpha value for each single point, RGBA = 4 components unsigned char rgba[4] = {0,0,0,0}; m_FiberColors = vtkSmartPointer::New(); m_FiberColors->Allocate(numOfPoints * 4); m_FiberColors->SetNumberOfComponents(4); m_FiberColors->SetName("FIBER_COLORS"); auto numOfFibers = m_FiberPolyData->GetNumberOfLines(); if (numOfFibers < 1) return; /* extract single fibers of fiberBundle */ vtkCellArray* fiberList = m_FiberPolyData->GetLines(); fiberList->InitTraversal(); for (int fi=0; fiGetNextCell(pointsPerFiber, idList); /* single fiber checkpoints: is number of points valid */ if (pointsPerFiber > 1) { /* operate on points of single fiber */ for (int i=0; i 0) { /* The color value of the current point is influenced by the previous point and next point. */ vnl_vector_fixed< double, 3 > currentPntvtk(extrPoints->GetPoint(idList[i])[0], extrPoints->GetPoint(idList[i])[1],extrPoints->GetPoint(idList[i])[2]); vnl_vector_fixed< double, 3 > nextPntvtk(extrPoints->GetPoint(idList[i+1])[0], extrPoints->GetPoint(idList[i+1])[1], extrPoints->GetPoint(idList[i+1])[2]); vnl_vector_fixed< double, 3 > prevPntvtk(extrPoints->GetPoint(idList[i-1])[0], extrPoints->GetPoint(idList[i-1])[1], extrPoints->GetPoint(idList[i-1])[2]); vnl_vector_fixed< double, 3 > diff1; diff1 = currentPntvtk - nextPntvtk; vnl_vector_fixed< double, 3 > diff2; diff2 = currentPntvtk - prevPntvtk; vnl_vector_fixed< double, 3 > diff; diff = (diff1 - diff2) / 2.0; diff.normalize(); rgba[0] = static_cast(255.0 * std::fabs(diff[0])); rgba[1] = static_cast(255.0 * std::fabs(diff[1])); rgba[2] = static_cast(255.0 * std::fabs(diff[2])); rgba[3] = static_cast(255.0); } else if (i==0) { /* First point has no previous point, therefore only diff1 is taken */ vnl_vector_fixed< double, 3 > currentPntvtk(extrPoints->GetPoint(idList[i])[0], extrPoints->GetPoint(idList[i])[1],extrPoints->GetPoint(idList[i])[2]); vnl_vector_fixed< double, 3 > nextPntvtk(extrPoints->GetPoint(idList[i+1])[0], extrPoints->GetPoint(idList[i+1])[1], extrPoints->GetPoint(idList[i+1])[2]); vnl_vector_fixed< double, 3 > diff1; diff1 = currentPntvtk - nextPntvtk; diff1.normalize(); rgba[0] = static_cast(255.0 * std::fabs(diff1[0])); rgba[1] = static_cast(255.0 * std::fabs(diff1[1])); rgba[2] = static_cast(255.0 * std::fabs(diff1[2])); rgba[3] = static_cast(255.0); } else if (i==pointsPerFiber-1) { /* Last point has no next point, therefore only diff2 is taken */ vnl_vector_fixed< double, 3 > currentPntvtk(extrPoints->GetPoint(idList[i])[0], extrPoints->GetPoint(idList[i])[1],extrPoints->GetPoint(idList[i])[2]); vnl_vector_fixed< double, 3 > prevPntvtk(extrPoints->GetPoint(idList[i-1])[0], extrPoints->GetPoint(idList[i-1])[1], extrPoints->GetPoint(idList[i-1])[2]); vnl_vector_fixed< double, 3 > diff2; diff2 = currentPntvtk - prevPntvtk; diff2.normalize(); rgba[0] = static_cast(255.0 * std::fabs(diff2[0])); rgba[1] = static_cast(255.0 * std::fabs(diff2[1])); rgba[2] = static_cast(255.0 * std::fabs(diff2[2])); rgba[3] = static_cast(255.0); } m_FiberColors->InsertTypedTuple(idList[i], rgba); } } else if (pointsPerFiber == 1) { /* a single point does not define a fiber (use vertex mechanisms instead */ continue; } else { MITK_DEBUG << "Fiber with 0 points detected... please check your tractography algorithm!" ; continue; } } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::ColorFibersByCurvature(bool opacity, bool weight_fibers, mitk::LookupTable::LookupTableType type) { double window = 5; //colors and alpha value for each single point, RGBA = 4 components unsigned char rgba[4] = {0,0,0,0}; m_FiberColors = vtkSmartPointer::New(); m_FiberColors->Allocate(m_FiberPolyData->GetNumberOfPoints() * 4); m_FiberColors->SetNumberOfComponents(4); m_FiberColors->SetName("FIBER_COLORS"); std::vector< double > values; double min = 1; double max = 0; MITK_INFO << "Coloring fibers by curvature"; boost::progress_display disp(static_cast(m_FiberPolyData->GetNumberOfCells())); for (int i=0; iGetNumberOfCells(); i++) { ++disp; vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); double mean_curv = 0; // calculate curvatures for (int j=0; j > vectors; vnl_vector_fixed< double, 3 > meanV; meanV.fill(0.0); while(dist1) { double p1[3]; points->GetPoint(c-1, p1); double p2[3]; points->GetPoint(c, p2); vnl_vector_fixed< double, 3 > v; v[0] = p2[0]-p1[0]; v[1] = p2[1]-p1[1]; v[2] = p2[2]-p1[2]; dist += v.magnitude(); v.normalize(); vectors.push_back(v); meanV += v; c--; } c = j; dist = 0; while(distGetPoint(c, p1); double p2[3]; points->GetPoint(c+1, p2); vnl_vector_fixed< double, 3 > v; v[0] = p2[0]-p1[0]; v[1] = p2[1]-p1[1]; v[2] = p2[2]-p1[2]; dist += v.magnitude(); v.normalize(); vectors.push_back(v); meanV += v; c++; } meanV.normalize(); double dev = 0; for (unsigned int c=0; c1.0) angle = 1.0; if (angle<-1.0) angle = -1.0; dev += acos(angle)*180/itk::Math::pi; } if (vectors.size()>0) dev /= vectors.size(); if (weight_fibers) mean_curv += dev; dev = 1.0-dev/180.0; values.push_back(dev); if (devmax) max = dev; } if (weight_fibers) this->SetFiberWeight(i, mean_curv/numPoints); } mitk::LookupTable::Pointer mitkLookup = mitk::LookupTable::New(); mitkLookup->SetType(type); if (type!=mitk::LookupTable::MULTILABEL) mitkLookup->GetVtkLookupTable()->SetTableRange(min, max); unsigned int count = 0; for (int i=0; iGetNumberOfCells(); i++) { vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); for (int j=0; jGetColor(dev, color); rgba[0] = static_cast(255.0 * color[0]); rgba[1] = static_cast(255.0 * color[1]); rgba[2] = static_cast(255.0 * color[2]); if (opacity) rgba[3] = static_cast(255.0f * dev/max); else rgba[3] = static_cast(255.0); m_FiberColors->InsertTypedTuple(cell->GetPointId(j), rgba); count++; } } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::SetFiberOpacity(vtkDoubleArray* FAValArray) { for(long i=0; iGetNumberOfTuples(); i++) { double faValue = FAValArray->GetValue(i); faValue = faValue * 255.0; m_FiberColors->SetComponent(i,3, static_cast(faValue) ); } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::ResetFiberOpacity() { for(long i=0; iGetNumberOfTuples(); i++) m_FiberColors->SetComponent(i,3, 255.0 ); m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::ColorFibersByScalarMap(mitk::Image::Pointer FAimage, bool opacity, bool weight_fibers, mitk::LookupTable::LookupTableType type, double max_cap, bool interpolate) { if (FAimage->GetPixelType().GetComponentTypeAsString()=="unsigned char" || FAimage->GetPixelType().GetComponentTypeAsString()=="char" || FAimage->GetPixelType().GetComponentTypeAsString()=="long" || FAimage->GetPixelType().GetComponentTypeAsString()=="unsigned long" || FAimage->GetPixelType().GetComponentTypeAsString()=="short" || FAimage->GetPixelType().GetComponentTypeAsString()=="unsigned short" || FAimage->GetPixelType().GetComponentTypeAsString()=="unsigned int" || FAimage->GetPixelType().GetComponentTypeAsString()=="int") { typedef itk::Image ImageType; ImageType::Pointer itkImage = ImageType::New(); CastToItkImage(FAimage, itkImage); ColorFibersByScalarMap(itkImage, opacity, weight_fibers, type, max_cap, interpolate ); } else { typedef itk::Image ImageType; ImageType::Pointer itkImage = ImageType::New(); CastToItkImage(FAimage, itkImage); ColorFibersByScalarMap(itkImage, opacity, weight_fibers, type, max_cap, interpolate ); } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } template void mitk::FiberBundle::ColorFibersByScalarMap(typename itk::Image::Pointer image, bool opacity, bool weight_fibers, mitk::LookupTable::LookupTableType type, double max_cap, bool interpolate) { m_FiberColors = vtkSmartPointer::New(); m_FiberColors->Allocate(m_FiberPolyData->GetNumberOfPoints() * 4); m_FiberColors->SetNumberOfComponents(4); m_FiberColors->SetName("FIBER_COLORS"); unsigned char rgba[4] = {0,0,0,0}; vtkPoints* pointSet = m_FiberPolyData->GetPoints(); if (type==mitk::LookupTable::MULTILABEL) interpolate = false; auto interpolator = itk::LinearInterpolateImageFunction< itk::Image, float >::New(); interpolator->SetInputImage(image); double min = 999999; double max = -999999; for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); double mean_val = 0; for (int j=0; jGetPoint(j, p); auto pixelValue = mitk::imv::GetImageValue(mitk::imv::GetItkPoint(p), interpolate, interpolator); if (pixelValue>max) max = pixelValue; if (pixelValueSetFiberWeight(i, mean_val/numPoints); } mitk::LookupTable::Pointer mitkLookup = mitk::LookupTable::New(); mitkLookup->SetType(type); if (type!=mitk::LookupTable::MULTILABEL) mitkLookup->GetVtkLookupTable()->SetTableRange(min, max*max_cap); for(long i=0; iGetNumberOfPoints(); ++i) { itk::Point px; px[0] = pointSet->GetPoint(i)[0]; px[1] = pointSet->GetPoint(i)[1]; px[2] = pointSet->GetPoint(i)[2]; auto pixelValue = mitk::imv::GetImageValue(px, interpolate, interpolator); double color[3]; mitkLookup->GetColor(pixelValue, color); rgba[0] = static_cast(255.0 * color[0]); rgba[1] = static_cast(255.0 * color[1]); rgba[2] = static_cast(255.0 * color[2]); if (opacity) rgba[3] = static_cast(255.0 * pixelValue); else rgba[3] = static_cast(255.0); m_FiberColors->InsertTypedTuple(i, rgba); } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::ColorFibersByFiberWeights(bool opacity, mitk::LookupTable::LookupTableType type) { m_FiberColors = vtkSmartPointer::New(); m_FiberColors->Allocate(m_FiberPolyData->GetNumberOfPoints() * 4); m_FiberColors->SetNumberOfComponents(4); m_FiberColors->SetName("FIBER_COLORS"); unsigned char rgba[4] = {0,0,0,0}; unsigned int counter = 0; float max = -999999; float min = 999999; for (unsigned int i=0; iGetFiberWeight(i); if (weight>max) max = weight; if (weightSetType(type); if (type!=mitk::LookupTable::MULTILABEL) mitkLookup->GetVtkLookupTable()->SetTableRange(min, max); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); auto weight = this->GetFiberWeight(i); double color[3]; mitkLookup->GetColor(weight, color); for (int j=0; j(255.0 * color[0]); rgba[1] = static_cast(255.0 * color[1]); rgba[2] = static_cast(255.0 * color[2]); if (opacity) rgba[3] = static_cast(255.0f * weight/max); else rgba[3] = static_cast(255.0); m_FiberColors->InsertTypedTuple(counter, rgba); counter++; } } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::SetSingleFiberColor(float r, float g, float b, unsigned int cellId, float alpha) { // if (m_FiberColors==nullptr) m_FiberColors = vtkSmartPointer::New(); m_FiberColors->Allocate(m_FiberPolyData->GetNumberOfPoints() * 4); m_FiberColors->SetNumberOfComponents(4); m_FiberColors->SetName("FIBER_COLORS"); // MITK_INFO << color->GetNumberOfTuples(); // MITK_INFO << m_FiberColors->GetNumberOfTuples(); unsigned char rgba[4] = {0,0,0,0}; unsigned int counter = 0; for (unsigned int i=0; iGetNumberOfCells(); i++) { vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); if (i==cellId) { - MITK_INFO << i; for (int j=0; j(r); rgba[1] = static_cast(b); rgba[2] = static_cast(g); rgba[3] = static_cast(alpha); // m_FiberColors->InsertTypedTuple(j, rgba); m_FiberColors->InsertTypedTuple(counter, rgba); counter++; } } else { for (int j=0; j(255); rgba[1] = static_cast(255); rgba[2] = static_cast(255); rgba[3] = static_cast(alpha); // m_FiberColors->InsertTypedTuple(j, rgba); m_FiberColors->InsertTypedTuple(counter, rgba); counter++; } } } // MITK_INFO << m_FiberColors->GetNumberOfTuples(); m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::SetFiberColors(float r, float g, float b, float alpha) { m_FiberColors = vtkSmartPointer::New(); m_FiberColors->Allocate(m_FiberPolyData->GetNumberOfPoints() * 4); m_FiberColors->SetNumberOfComponents(4); m_FiberColors->SetName("FIBER_COLORS"); unsigned char rgba[4] = {0,0,0,0}; for(long i=0; iGetNumberOfPoints(); ++i) { rgba[0] = static_cast(r); rgba[1] = static_cast(g); rgba[2] = static_cast(b); rgba[3] = static_cast(alpha); m_FiberColors->InsertTypedTuple(i, rgba); } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } void mitk::FiberBundle::GenerateFiberIds() { if (m_FiberPolyData == nullptr) return; vtkSmartPointer idFiberFilter = vtkSmartPointer::New(); idFiberFilter->SetInputData(m_FiberPolyData); idFiberFilter->CellIdsOn(); // idFiberFilter->PointIdsOn(); // point id's are not needed idFiberFilter->SetCellIdsArrayName(FIBER_ID_ARRAY); idFiberFilter->FieldDataOn(); idFiberFilter->Update(); m_FiberIdDataSet = idFiberFilter->GetOutput(); } float mitk::FiberBundle::GetNumEpFractionInMask(ItkUcharImgType* mask, bool different_label) { vtkSmartPointer PolyData = m_FiberPolyData; MITK_INFO << "Calculating EP-Fraction"; boost::progress_display disp(m_NumFibers); unsigned int in_mask = 0; for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); itk::Point startVertex =mitk::imv::GetItkPoint(points->GetPoint(0)); itk::Index<3> startIndex; mask->TransformPhysicalPointToIndex(startVertex, startIndex); itk::Point endVertex =mitk::imv::GetItkPoint(points->GetPoint(numPoints-1)); itk::Index<3> endIndex; mask->TransformPhysicalPointToIndex(endVertex, endIndex); if (mask->GetLargestPossibleRegion().IsInside(startIndex) && mask->GetLargestPossibleRegion().IsInside(endIndex)) { float v1 = mask->GetPixel(startIndex); if (v1 < 0.5f) continue; float v2 = mask->GetPixel(startIndex); if (v2 < 0.5f) continue; if (!different_label) ++in_mask; else if (fabs(v1-v2)>0.00001f) ++in_mask; } } return float(in_mask)/m_NumFibers; } std::tuple mitk::FiberBundle::GetDirectionalOverlap(ItkUcharImgType* mask, mitk::PeakImage::ItkPeakImageType* peak_image) { vtkSmartPointer PolyData = m_FiberPolyData; MITK_INFO << "Calculating overlap"; auto spacing = mask->GetSpacing(); boost::progress_display disp(m_NumFibers); double length_sum = 0; double in_mask_length = 0; double aligned_length = 0; for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); for (int j=0; j startVertex =mitk::imv::GetItkPoint(points->GetPoint(j)); itk::Index<3> startIndex; itk::ContinuousIndex startIndexCont; mask->TransformPhysicalPointToIndex(startVertex, startIndex); mask->TransformPhysicalPointToContinuousIndex(startVertex, startIndexCont); itk::Point endVertex =mitk::imv::GetItkPoint(points->GetPoint(j + 1)); itk::Index<3> endIndex; itk::ContinuousIndex endIndexCont; mask->TransformPhysicalPointToIndex(endVertex, endIndex); mask->TransformPhysicalPointToContinuousIndex(endVertex, endIndexCont); vnl_vector_fixed< float, 3 > fdir; fdir[0] = endVertex[0] - startVertex[0]; fdir[1] = endVertex[1] - startVertex[1]; fdir[2] = endVertex[2] - startVertex[2]; fdir.normalize(); std::vector< std::pair< itk::Index<3>, double > > segments = mitk::imv::IntersectImage(spacing, startIndex, endIndex, startIndexCont, endIndexCont); for (std::pair< itk::Index<3>, double > segment : segments) { if ( mask->GetLargestPossibleRegion().IsInside(segment.first) && mask->GetPixel(segment.first) > 0 ) { in_mask_length += segment.second; mitk::PeakImage::ItkPeakImageType::IndexType idx4; idx4[0] = segment.first[0]; idx4[1] = segment.first[1]; idx4[2] = segment.first[2]; vnl_vector_fixed< float, 3 > peak; idx4[3] = 0; peak[0] = peak_image->GetPixel(idx4); idx4[3] = 1; peak[1] = peak_image->GetPixel(idx4); idx4[3] = 2; peak[2] = peak_image->GetPixel(idx4); if (std::isnan(peak[0]) || std::isnan(peak[1]) || std::isnan(peak[2]) || peak.magnitude()<0.0001f) continue; peak.normalize(); double f = 1.0 - std::acos(std::fabs(static_cast(dot_product(fdir, peak)))) * 2.0/itk::Math::pi; aligned_length += segment.second * f; } length_sum += segment.second; } } } if (length_sum<=0.0001) { MITK_INFO << "Fiber length sum is zero!"; return std::make_tuple(0,0); } return std::make_tuple(aligned_length/length_sum, in_mask_length/length_sum); } float mitk::FiberBundle::GetOverlap(ItkUcharImgType* mask) { vtkSmartPointer PolyData = m_FiberPolyData; MITK_INFO << "Calculating overlap"; auto spacing = mask->GetSpacing(); boost::progress_display disp(m_NumFibers); double length_sum = 0; double in_mask_length = 0; for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); for (int j=0; j startVertex =mitk::imv::GetItkPoint(points->GetPoint(j)); itk::Index<3> startIndex; itk::ContinuousIndex startIndexCont; mask->TransformPhysicalPointToIndex(startVertex, startIndex); mask->TransformPhysicalPointToContinuousIndex(startVertex, startIndexCont); itk::Point endVertex =mitk::imv::GetItkPoint(points->GetPoint(j + 1)); itk::Index<3> endIndex; itk::ContinuousIndex endIndexCont; mask->TransformPhysicalPointToIndex(endVertex, endIndex); mask->TransformPhysicalPointToContinuousIndex(endVertex, endIndexCont); std::vector< std::pair< itk::Index<3>, double > > segments = mitk::imv::IntersectImage(spacing, startIndex, endIndex, startIndexCont, endIndexCont); for (std::pair< itk::Index<3>, double > segment : segments) { if ( mask->GetLargestPossibleRegion().IsInside(segment.first) && mask->GetPixel(segment.first) > 0 ) in_mask_length += segment.second; length_sum += segment.second; } } } if (length_sum<=0.000001) { MITK_INFO << "Fiber length sum is zero!"; return 0; } return static_cast(in_mask_length/length_sum); } mitk::FiberBundle::Pointer mitk::FiberBundle::RemoveFibersOutside(ItkUcharImgType* mask, bool invert) { vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); std::vector< float > fib_weights; MITK_INFO << "Cutting fibers"; boost::progress_display disp(m_NumFibers); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); int newNumPoints = 0; if (numPoints>1) { for (int j=0; j itkP =mitk::imv::GetItkPoint(points->GetPoint(j)); itk::Index<3> idx; mask->TransformPhysicalPointToIndex(itkP, idx); bool inside = false; if ( mask->GetLargestPossibleRegion().IsInside(idx) && mask->GetPixel(idx)!=0 ) inside = true; if (inside && !invert) { vtkIdType id = vtkNewPoints->InsertNextPoint(itkP.GetDataPointer()); container->GetPointIds()->InsertNextId(id); newNumPoints++; } else if ( !inside && invert ) { vtkIdType id = vtkNewPoints->InsertNextPoint(itkP.GetDataPointer()); container->GetPointIds()->InsertNextId(id); newNumPoints++; } else if (newNumPoints>1) { fib_weights.push_back(this->GetFiberWeight(i)); vtkNewCells->InsertNextCell(container); newNumPoints = 0; container = vtkSmartPointer::New(); } else { newNumPoints = 0; container = vtkSmartPointer::New(); } } if (newNumPoints>1) { fib_weights.push_back(this->GetFiberWeight(i)); vtkNewCells->InsertNextCell(container); } } } vtkSmartPointer newFiberWeights = vtkSmartPointer::New(); newFiberWeights->SetName("FIBER_WEIGHTS"); newFiberWeights->SetNumberOfValues(static_cast(fib_weights.size())); if (vtkNewCells->GetNumberOfCells()<=0) return nullptr; for (unsigned int i=0; iGetNumberOfValues(); i++) newFiberWeights->SetValue(i, fib_weights.at(i)); // vtkSmartPointer newFiberColors = vtkSmartPointer::New(); // newFiberColors->Allocate(m_FiberPolyData->GetNumberOfPoints() * 4); // newFiberColors->SetNumberOfComponents(4); // newFiberColors->SetName("FIBER_COLORS"); // unsigned char rgba[4] = {0,0,0,0}; // for(long i=0; iGetNumberOfPoints(); ++i) // { // rgba[0] = (unsigned char) r; // rgba[1] = (unsigned char) g; // rgba[2] = (unsigned char) b; // rgba[3] = (unsigned char) alpha; // m_FiberColors->InsertTypedTuple(i, rgba); // } vtkSmartPointer newPolyData = vtkSmartPointer::New(); newPolyData->SetPoints(vtkNewPoints); newPolyData->SetLines(vtkNewCells); mitk::FiberBundle::Pointer newFib = mitk::FiberBundle::New(newPolyData); newFib->SetFiberWeights(newFiberWeights); // newFib->Compress(0.1); newFib->SetTrackVisHeader(this->GetTrackVisHeader()); return newFib; } mitk::FiberBundle::Pointer mitk::FiberBundle::ExtractFiberSubset(DataNode* roi, DataStorage* storage) { if (roi==nullptr || !(dynamic_cast(roi->GetData()) || dynamic_cast(roi->GetData())) ) return nullptr; std::vector tmp = ExtractFiberIdSubset(roi, storage); if (tmp.size()<=0) return mitk::FiberBundle::New(); vtkSmartPointer weights = vtkSmartPointer::New(); vtkSmartPointer pTmp = GeneratePolyDataByIds(tmp, weights); mitk::FiberBundle::Pointer fib = mitk::FiberBundle::New(pTmp); fib->SetFiberWeights(weights); fib->SetTrackVisHeader(this->GetTrackVisHeader()); return fib; } std::vector mitk::FiberBundle::ExtractFiberIdSubset(DataNode *roi, DataStorage* storage) { std::vector result; if (roi==nullptr || roi->GetData()==nullptr) return result; mitk::PlanarFigureComposite::Pointer pfc = dynamic_cast(roi->GetData()); if (!pfc.IsNull()) // handle composite { DataStorage::SetOfObjects::ConstPointer children = storage->GetDerivations(roi); if (children->size()==0) return result; switch (pfc->getOperationType()) { case 0: // AND { MITK_INFO << "AND"; result = this->ExtractFiberIdSubset(children->ElementAt(0), storage); std::vector::iterator it; for (unsigned int i=1; iSize(); ++i) { std::vector inRoi = this->ExtractFiberIdSubset(children->ElementAt(i), storage); std::vector rest(std::min(result.size(),inRoi.size())); it = std::set_intersection(result.begin(), result.end(), inRoi.begin(), inRoi.end(), rest.begin() ); rest.resize( static_cast(it - rest.begin()) ); result = rest; } break; } case 1: // OR { MITK_INFO << "OR"; result = ExtractFiberIdSubset(children->ElementAt(0), storage); std::vector::iterator it; for (unsigned int i=1; iSize(); ++i) { it = result.end(); std::vector inRoi = ExtractFiberIdSubset(children->ElementAt(i), storage); result.insert(it, inRoi.begin(), inRoi.end()); } // remove duplicates sort(result.begin(), result.end()); it = unique(result.begin(), result.end()); result.resize( static_cast(it - result.begin()) ); break; } case 2: // NOT { MITK_INFO << "NOT"; for(unsigned int i=0; iGetNumFibers(); i++) result.push_back(i); std::vector::iterator it; for (unsigned int i=0; iSize(); ++i) { std::vector inRoi = ExtractFiberIdSubset(children->ElementAt(i), storage); std::vector rest(result.size()-inRoi.size()); it = std::set_difference(result.begin(), result.end(), inRoi.begin(), inRoi.end(), rest.begin() ); rest.resize( static_cast(it - rest.begin()) ); result = rest; } break; } } } else if ( dynamic_cast(roi->GetData()) ) // actual extraction { if ( dynamic_cast(roi->GetData()) ) { mitk::PlanarFigure::Pointer planarPoly = dynamic_cast(roi->GetData()); //create vtkPolygon using controlpoints from planarFigure polygon vtkSmartPointer polygonVtk = vtkSmartPointer::New(); for (unsigned int i=0; iGetNumberOfControlPoints(); ++i) { itk::Point p = planarPoly->GetWorldControlPoint(i); vtkIdType id = polygonVtk->GetPoints()->InsertNextPoint(p[0], p[1], p[2] ); polygonVtk->GetPointIds()->InsertNextId(id); } MITK_INFO << "Extracting with polygon"; boost::progress_display disp(m_NumFibers); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); for (int j=0; jGetPoint(j, p1); double p2[3] = {0,0,0}; points->GetPoint(j+1, p2); double tolerance = 0.001; // Outputs double t = 0; // Parametric coordinate of intersection (0 (corresponding to p1) to 1 (corresponding to p2)) double x[3] = {0,0,0}; // The coordinate of the intersection double pcoords[3] = {0,0,0}; int subId = 0; int iD = polygonVtk->IntersectWithLine(p1, p2, tolerance, t, x, pcoords, subId); if (iD!=0) { result.push_back(i); break; } } } } else if ( dynamic_cast(roi->GetData()) ) { mitk::PlanarFigure::Pointer planarFigure = dynamic_cast(roi->GetData()); Vector3D planeNormal = planarFigure->GetPlaneGeometry()->GetNormal(); planeNormal.Normalize(); //calculate circle radius mitk::Point3D V1w = planarFigure->GetWorldControlPoint(0); //centerPoint mitk::Point3D V2w = planarFigure->GetWorldControlPoint(1); //radiusPoint double radius = V1w.EuclideanDistanceTo(V2w); radius *= radius; MITK_INFO << "Extracting with circle"; boost::progress_display disp(m_NumFibers); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); for (int j=0; jGetPoint(j, p1); double p2[3] = {0,0,0}; points->GetPoint(j+1, p2); // Outputs double t = 0; // Parametric coordinate of intersection (0 (corresponding to p1) to 1 (corresponding to p2)) double x[3] = {0,0,0}; // The coordinate of the intersection int iD = vtkPlane::IntersectWithLine(p1,p2,planeNormal.GetDataPointer(),V1w.GetDataPointer(),t,x); if (iD!=0) { double dist = (x[0]-V1w[0])*(x[0]-V1w[0])+(x[1]-V1w[1])*(x[1]-V1w[1])+(x[2]-V1w[2])*(x[2]-V1w[2]); if( dist <= radius) { result.push_back(i); break; } } } } } return result; } return result; } void mitk::FiberBundle::UpdateFiberGeometry() { vtkSmartPointer cleaner = vtkSmartPointer::New(); cleaner->SetInputData(m_FiberPolyData); cleaner->PointMergingOff(); cleaner->Update(); m_FiberPolyData = cleaner->GetOutput(); m_FiberLengths.clear(); m_MeanFiberLength = 0; m_MedianFiberLength = 0; m_LengthStDev = 0; m_NumFibers = static_cast(m_FiberPolyData->GetNumberOfCells()); if (m_FiberColors==nullptr || m_FiberColors->GetNumberOfTuples()!=m_FiberPolyData->GetNumberOfPoints()) this->ColorFibersByOrientation(); if (m_FiberWeights->GetNumberOfValues()!=m_NumFibers) { m_FiberWeights = vtkSmartPointer::New(); m_FiberWeights->SetName("FIBER_WEIGHTS"); m_FiberWeights->SetNumberOfValues(m_NumFibers); this->SetFiberWeights(1); } if (m_NumFibers<=0) // no fibers present; apply default geometry { m_MinFiberLength = 0; m_MaxFiberLength = 0; mitk::Geometry3D::Pointer geometry = mitk::Geometry3D::New(); geometry->SetImageGeometry(false); float b[] = {0, 1, 0, 1, 0, 1}; geometry->SetFloatBounds(b); SetGeometry(geometry); return; } double b[6]; m_FiberPolyData->GetBounds(b); // calculate statistics for (int i=0; iGetNumberOfCells(); i++) { vtkCell* cell = m_FiberPolyData->GetCell(i); auto p = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); float length = 0; for (int j=0; jGetPoint(j, p1); double p2[3]; points->GetPoint(j+1, p2); double dist = std::sqrt((p1[0]-p2[0])*(p1[0]-p2[0])+(p1[1]-p2[1])*(p1[1]-p2[1])+(p1[2]-p2[2])*(p1[2]-p2[2])); length += static_cast(dist); } m_FiberLengths.push_back(length); m_MeanFiberLength += length; if (i==0) { m_MinFiberLength = length; m_MaxFiberLength = length; } else { if (lengthm_MaxFiberLength) m_MaxFiberLength = length; } } m_MeanFiberLength /= m_NumFibers; std::vector< float > sortedLengths = m_FiberLengths; std::sort(sortedLengths.begin(), sortedLengths.end()); for (unsigned int i=0; i1) m_LengthStDev /= (m_NumFibers-1); else m_LengthStDev = 0; m_LengthStDev = std::sqrt(m_LengthStDev); m_MedianFiberLength = sortedLengths.at(m_NumFibers/2); mitk::Geometry3D::Pointer geometry = mitk::Geometry3D::New(); geometry->SetFloatBounds(b); this->SetGeometry(geometry); GetTrackVisHeader(); m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } float mitk::FiberBundle::GetFiberWeight(unsigned int fiber) const { return m_FiberWeights->GetValue(fiber); } void mitk::FiberBundle::SetFiberWeights(float newWeight) { for (int i=0; iGetNumberOfValues(); i++) m_FiberWeights->SetValue(i, newWeight); } void mitk::FiberBundle::SetFiberWeights(vtkSmartPointer weights) { if (m_NumFibers!=weights->GetNumberOfValues()) { MITK_INFO << "Weights array not equal to number of fibers! " << weights->GetNumberOfValues() << " vs " << m_NumFibers; return; } for (int i=0; iGetNumberOfValues(); i++) m_FiberWeights->SetValue(i, weights->GetValue(i)); m_FiberWeights->SetName("FIBER_WEIGHTS"); } void mitk::FiberBundle::SetFiberWeight(unsigned int fiber, float weight) { m_FiberWeights->SetValue(fiber, weight); } void mitk::FiberBundle::SetFiberColors(vtkSmartPointer fiberColors) { for(long i=0; iGetNumberOfPoints(); ++i) { unsigned char source[4] = {0,0,0,0}; fiberColors->GetTypedTuple(i, source); unsigned char target[4] = {0,0,0,0}; target[0] = source[0]; target[1] = source[1]; target[2] = source[2]; target[3] = source[3]; m_FiberColors->InsertTypedTuple(i, target); } m_UpdateTime3D.Modified(); m_UpdateTime2D.Modified(); } itk::Matrix< double, 3, 3 > mitk::FiberBundle::TransformMatrix(itk::Matrix< double, 3, 3 > m, double rx, double ry, double rz) { rx = rx*itk::Math::pi/180; ry = ry*itk::Math::pi/180; rz = rz*itk::Math::pi/180; itk::Matrix< double, 3, 3 > rotX; rotX.SetIdentity(); rotX[1][1] = cos(rx); rotX[2][2] = rotX[1][1]; rotX[1][2] = -sin(rx); rotX[2][1] = -rotX[1][2]; itk::Matrix< double, 3, 3 > rotY; rotY.SetIdentity(); rotY[0][0] = cos(ry); rotY[2][2] = rotY[0][0]; rotY[0][2] = sin(ry); rotY[2][0] = -rotY[0][2]; itk::Matrix< double, 3, 3 > rotZ; rotZ.SetIdentity(); rotZ[0][0] = cos(rz); rotZ[1][1] = rotZ[0][0]; rotZ[0][1] = -sin(rz); rotZ[1][0] = -rotZ[0][1]; itk::Matrix< double, 3, 3 > rot = rotZ*rotY*rotX; m = rot*m; return m; } void mitk::FiberBundle::TransformFibers(itk::ScalableAffineTransform< mitk::ScalarType >::Pointer transform) { vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; j p =mitk::imv::GetItkPoint(points->GetPoint(j)); p = transform->TransformPoint(p); vtkIdType id = vtkNewPoints->InsertNextPoint(p.GetDataPointer()); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); } m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } void mitk::FiberBundle::TransformFibers(double rx, double ry, double rz, double tx, double ty, double tz) { vnl_matrix_fixed< double, 3, 3 > rot = mitk::imv::GetRotationMatrixVnl(rx, ry, rz); mitk::BaseGeometry::Pointer geom = this->GetGeometry(); mitk::Point3D center = geom->GetCenter(); vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j); vnl_vector_fixed< double, 3 > dir; dir[0] = p[0]-center[0]; dir[1] = p[1]-center[1]; dir[2] = p[2]-center[2]; dir = rot*dir; dir[0] += center[0]+tx; dir[1] += center[1]+ty; dir[2] += center[2]+tz; vtkIdType id = vtkNewPoints->InsertNextPoint(dir.data_block()); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); } m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } void mitk::FiberBundle::RotateAroundAxis(double x, double y, double z) { x = x*itk::Math::pi/180; y = y*itk::Math::pi/180; z = z*itk::Math::pi/180; vnl_matrix_fixed< double, 3, 3 > rotX; rotX.set_identity(); rotX[1][1] = cos(x); rotX[2][2] = rotX[1][1]; rotX[1][2] = -sin(x); rotX[2][1] = -rotX[1][2]; vnl_matrix_fixed< double, 3, 3 > rotY; rotY.set_identity(); rotY[0][0] = cos(y); rotY[2][2] = rotY[0][0]; rotY[0][2] = sin(y); rotY[2][0] = -rotY[0][2]; vnl_matrix_fixed< double, 3, 3 > rotZ; rotZ.set_identity(); rotZ[0][0] = cos(z); rotZ[1][1] = rotZ[0][0]; rotZ[0][1] = -sin(z); rotZ[1][0] = -rotZ[0][1]; mitk::BaseGeometry::Pointer geom = this->GetGeometry(); mitk::Point3D center = geom->GetCenter(); vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j); vnl_vector_fixed< double, 3 > dir; dir[0] = p[0]-center[0]; dir[1] = p[1]-center[1]; dir[2] = p[2]-center[2]; dir = rotZ*rotY*rotX*dir; dir[0] += center[0]; dir[1] += center[1]; dir[2] += center[2]; vtkIdType id = vtkNewPoints->InsertNextPoint(dir.data_block()); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); } m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } void mitk::FiberBundle::ScaleFibers(double x, double y, double z, bool subtractCenter) { MITK_INFO << "Scaling fibers"; boost::progress_display disp(m_NumFibers); mitk::BaseGeometry* geom = this->GetGeometry(); mitk::Point3D c = geom->GetCenter(); vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j); if (subtractCenter) { p[0] -= c[0]; p[1] -= c[1]; p[2] -= c[2]; } p[0] *= x; p[1] *= y; p[2] *= z; if (subtractCenter) { p[0] += c[0]; p[1] += c[1]; p[2] += c[2]; } vtkIdType id = vtkNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); } m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } void mitk::FiberBundle::TranslateFibers(double x, double y, double z) { vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j); p[0] += x; p[1] += y; p[2] += z; vtkIdType id = vtkNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); } m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } void mitk::FiberBundle::MirrorFibers(unsigned int axis) { if (axis>2) return; MITK_INFO << "Mirroring fibers"; boost::progress_display disp(m_NumFibers); vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j); p[axis] = -p[axis]; vtkIdType id = vtkNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); } m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } void mitk::FiberBundle::RemoveDir(vnl_vector_fixed dir, double threshold) { dir.normalize(); vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); boost::progress_display disp(static_cast(m_FiberPolyData->GetNumberOfCells())); for (int i=0; iGetNumberOfCells(); i++) { ++disp ; vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); // calculate curvatures vtkSmartPointer container = vtkSmartPointer::New(); bool discard = false; for (int j=0; jGetPoint(j, p1); double p2[3]; points->GetPoint(j+1, p2); vnl_vector_fixed< double, 3 > v1; v1[0] = p2[0]-p1[0]; v1[1] = p2[1]-p1[1]; v1[2] = p2[2]-p1[2]; if (v1.magnitude()>0.001) { v1.normalize(); if (fabs(dot_product(v1,dir))>threshold) { discard = true; break; } } } if (!discard) { for (int j=0; jGetPoint(j, p1); vtkIdType id = vtkNewPoints->InsertNextPoint(p1); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); } } m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); // UpdateColorCoding(); // UpdateFiberGeometry(); } bool mitk::FiberBundle::ApplyCurvatureThreshold(float minRadius, bool deleteFibers) { if (minRadius<0) return true; vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); MITK_INFO << "Applying curvature threshold"; boost::progress_display disp(static_cast(m_FiberPolyData->GetNumberOfCells())); for (int i=0; iGetNumberOfCells(); i++) { ++disp ; vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); // calculate curvatures vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j, p1); double p2[3]; points->GetPoint(j+1, p2); double p3[3]; points->GetPoint(j+2, p3); vnl_vector_fixed< float, 3 > v1, v2, v3; v1[0] = static_cast(p2[0]-p1[0]); v1[1] = static_cast(p2[1]-p1[1]); v1[2] = static_cast(p2[2]-p1[2]); v2[0] = static_cast(p3[0]-p2[0]); v2[1] = static_cast(p3[1]-p2[1]); v2[2] = static_cast(p3[2]-p2[2]); v3[0] = static_cast(p1[0]-p3[0]); v3[1] = static_cast(p1[1]-p3[1]); v3[2] = static_cast(p1[2]-p3[2]); float a = v1.magnitude(); float b = v2.magnitude(); float c = v3.magnitude(); float r = a*b*c/std::sqrt((a+b+c)*(a+b-c)*(b+c-a)*(a-b+c)); // radius of triangle via Heron's formula (area of triangle) vtkIdType id = vtkNewPoints->InsertNextPoint(p1); container->GetPointIds()->InsertNextId(id); if (deleteFibers && rInsertNextCell(container); container = vtkSmartPointer::New(); } else if (j==numPoints-3) { id = vtkNewPoints->InsertNextPoint(p2); container->GetPointIds()->InsertNextId(id); id = vtkNewPoints->InsertNextPoint(p3); container->GetPointIds()->InsertNextId(id); vtkNewCells->InsertNextCell(container); } } } if (vtkNewCells->GetNumberOfCells()<=0) return false; m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); return true; } bool mitk::FiberBundle::RemoveShortFibers(float lengthInMM) { MITK_INFO << "Removing short fibers"; if (lengthInMM<=0 || lengthInMMm_MaxFiberLength) // can't remove all fibers { MITK_WARN << "Process aborted. No fibers would be left!"; return false; } vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); float min = m_MaxFiberLength; boost::progress_display disp(m_NumFibers); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); if (m_FiberLengths.at(i)>=lengthInMM) { vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j); vtkIdType id = vtkNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); if (m_FiberLengths.at(i)GetNumberOfCells()<=0) return false; m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); return true; } bool mitk::FiberBundle::RemoveLongFibers(float lengthInMM) { if (lengthInMM<=0 || lengthInMM>m_MaxFiberLength) return true; if (lengthInMM vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); MITK_INFO << "Removing long fibers"; boost::progress_display disp(m_NumFibers); for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); if (m_FiberLengths.at(i)<=lengthInMM) { vtkSmartPointer container = vtkSmartPointer::New(); for (int j=0; jGetPoint(j); vtkIdType id = vtkNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vtkNewCells->InsertNextCell(container); } } if (vtkNewCells->GetNumberOfCells()<=0) return false; m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); return true; } void mitk::FiberBundle::ResampleSpline(float pointDistance, double tension, double continuity, double bias ) { if (pointDistance<=0) return; vtkSmartPointer vtkSmoothPoints = vtkSmartPointer::New(); //in smoothpoints the interpolated points representing a fiber are stored. //in vtkcells all polylines are stored, actually all id's of them are stored vtkSmartPointer vtkSmoothCells = vtkSmartPointer::New(); //cellcontainer for smoothed lines MITK_INFO << "Smoothing fibers"; vtkSmartPointer newFiberWeights = vtkSmartPointer::New(); newFiberWeights->SetName("FIBER_WEIGHTS"); newFiberWeights->SetNumberOfValues(m_NumFibers); std::vector< vtkSmartPointer > resampled_streamlines; resampled_streamlines.resize(m_NumFibers); boost::progress_display disp(m_NumFibers); #pragma omp parallel for for (int i=0; i(m_NumFibers); i++) { vtkSmartPointer newPoints = vtkSmartPointer::New(); float length = 0; #pragma omp critical { length = m_FiberLengths.at(static_cast(i)); ++disp; vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); for (int j=0; jInsertNextPoint(points->GetPoint(j)); } int sampling = static_cast(std::ceil(length/pointDistance)); vtkSmartPointer xSpline = vtkSmartPointer::New(); vtkSmartPointer ySpline = vtkSmartPointer::New(); vtkSmartPointer zSpline = vtkSmartPointer::New(); xSpline->SetDefaultBias(bias); xSpline->SetDefaultTension(tension); xSpline->SetDefaultContinuity(continuity); ySpline->SetDefaultBias(bias); ySpline->SetDefaultTension(tension); ySpline->SetDefaultContinuity(continuity); zSpline->SetDefaultBias(bias); zSpline->SetDefaultTension(tension); zSpline->SetDefaultContinuity(continuity); vtkSmartPointer spline = vtkSmartPointer::New(); spline->SetXSpline(xSpline); spline->SetYSpline(ySpline); spline->SetZSpline(zSpline); spline->SetPoints(newPoints); vtkSmartPointer functionSource = vtkSmartPointer::New(); functionSource->SetParametricFunction(spline); functionSource->SetUResolution(sampling); functionSource->SetVResolution(sampling); functionSource->SetWResolution(sampling); functionSource->Update(); vtkPolyData* outputFunction = functionSource->GetOutput(); vtkPoints* tmpSmoothPnts = outputFunction->GetPoints(); //smoothPoints of current fiber vtkSmartPointer smoothLine = vtkSmartPointer::New(); #pragma omp critical { for (int j=0; jGetNumberOfPoints(); j++) { vtkIdType id = vtkSmoothPoints->InsertNextPoint(tmpSmoothPnts->GetPoint(j)); smoothLine->GetPointIds()->InsertNextId(id); } resampled_streamlines[static_cast(i)] = smoothLine; } } for (auto container : resampled_streamlines) { vtkSmoothCells->InsertNextCell(container); } m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkSmoothPoints); m_FiberPolyData->SetLines(vtkSmoothCells); this->SetFiberPolyData(m_FiberPolyData, true); } void mitk::FiberBundle::ResampleSpline(float pointDistance) { ResampleSpline(pointDistance, 0, 0, 0 ); } unsigned int mitk::FiberBundle::GetNumberOfPoints() const { unsigned int points = 0; for (int i=0; iGetNumberOfCells(); i++) { vtkCell* cell = m_FiberPolyData->GetCell(i); points += cell->GetNumberOfPoints(); } return points; } void mitk::FiberBundle::Compress(float error) { vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); MITK_INFO << "Compressing fibers with max. error " << error << "mm"; unsigned int numRemovedPoints = 0; boost::progress_display disp(static_cast(m_FiberPolyData->GetNumberOfCells())); vtkSmartPointer newFiberWeights = vtkSmartPointer::New(); newFiberWeights->SetName("FIBER_WEIGHTS"); newFiberWeights->SetNumberOfValues(m_NumFibers); #pragma omp parallel for for (int i=0; i(m_FiberPolyData->GetNumberOfCells()); i++) { std::vector< vnl_vector_fixed< double, 3 > > vertices; float weight = 1; #pragma omp critical { ++disp; weight = m_FiberWeights->GetValue(i); vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); for (int j=0; jGetPoint(j, cand); vnl_vector_fixed< double, 3 > candV; candV[0]=cand[0]; candV[1]=cand[1]; candV[2]=cand[2]; vertices.push_back(candV); } } // calculate curvatures auto numPoints = vertices.size(); std::vector< int > removedPoints; removedPoints.resize(numPoints, 0); removedPoints[0]=-1; removedPoints[numPoints-1]=-1; vtkSmartPointer container = vtkSmartPointer::New(); unsigned int remCounter = 0; bool pointFound = true; while (pointFound) { pointFound = false; double minError = static_cast(error); unsigned int removeIndex = 0; for (unsigned int j=0; j candV = vertices.at(j); int validP = -1; vnl_vector_fixed< double, 3 > pred; for (int k=static_cast(j)-1; k>=0; k--) if (removedPoints[static_cast(k)]<=0) { pred = vertices.at(static_cast(k)); validP = k; break; } int validS = -1; vnl_vector_fixed< double, 3 > succ; for (unsigned int k=j+1; k(k); break; } if (validP>=0 && validS>=0) { double a = (candV-pred).magnitude(); double b = (candV-succ).magnitude(); double c = (pred-succ).magnitude(); double s=0.5*(a+b+c); double hc=(2.0/c)*sqrt(fabs(s*(s-a)*(s-b)*(s-c))); if (hcInsertNextPoint(vertices.at(j).data_block()); container->GetPointIds()->InsertNextId(id); } } } #pragma omp critical { newFiberWeights->SetValue(vtkNewCells->GetNumberOfCells(), weight); numRemovedPoints += remCounter; vtkNewCells->InsertNextCell(container); } } if (vtkNewCells->GetNumberOfCells()>0) { MITK_INFO << "Removed points: " << numRemovedPoints; SetFiberWeights(newFiberWeights); m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } } void mitk::FiberBundle::ResampleToNumPoints(unsigned int targetPoints) { if (targetPoints<2) mitkThrow() << "Minimum two points required for resampling!"; MITK_INFO << "Resampling fibers (number of points " << targetPoints << ")"; bool unequal_fibs = true; while (unequal_fibs) { vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); vtkSmartPointer newFiberWeights = vtkSmartPointer::New(); newFiberWeights->SetName("FIBER_WEIGHTS"); newFiberWeights->SetNumberOfValues(m_NumFibers); unequal_fibs = false; MITK_INFO << "Start"; // cv::parallel_for_(cv::Range(0, m_FiberPolyData->GetNumberOfCells()), [&](const cv::Range &range) // { // for (int i = range.start; i < range.end; i++) for (unsigned int i=0; iGetNumberOfCells(); i++) { std::vector< vnl_vector_fixed< double, 3 > > vertices; float weight = 1; double seg_len = 0; { weight = m_FiberWeights->GetValue(i); vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); if (numPoints!=targetPoints) seg_len = static_cast(this->GetFiberLength(i)/(targetPoints-1)); vtkPoints* points = cell->GetPoints(); for (int j=0; jGetPoint(j, cand); vnl_vector_fixed< double, 3 > candV; candV[0]=cand[0]; candV[1]=cand[1]; candV[2]=cand[2]; vertices.push_back(candV); } } vtkSmartPointer container = vtkSmartPointer::New(); vnl_vector_fixed< double, 3 > lastV = vertices.at(0); { vtkIdType id = vtkNewPoints->InsertNextPoint(lastV.data_block()); container->GetPointIds()->InsertNextId(id); } for (unsigned int j=1; j vec = vertices.at(j) - lastV; double new_dist = vec.magnitude(); if (new_dist >= seg_len && seg_len>0) { vnl_vector_fixed< double, 3 > newV = lastV; if ( new_dist-seg_len <= mitk::eps ) { vec.normalize(); newV += vec * seg_len; } else { // intersection between sphere (radius 'pointDistance', center 'lastV') and line (direction 'd' and point 'p') vnl_vector_fixed< double, 3 > p = vertices.at(j-1); vnl_vector_fixed< double, 3 > d = vertices.at(j) - p; double a = d[0]*d[0] + d[1]*d[1] + d[2]*d[2]; double b = 2 * (d[0] * (p[0] - lastV[0]) + d[1] * (p[1] - lastV[1]) + d[2] * (p[2] - lastV[2])); double c = (p[0] - lastV[0])*(p[0] - lastV[0]) + (p[1] - lastV[1])*(p[1] - lastV[1]) + (p[2] - lastV[2])*(p[2] - lastV[2]) - seg_len*seg_len; double v1 =(-b + std::sqrt(b*b-4*a*c))/(2*a); double v2 =(-b - std::sqrt(b*b-4*a*c))/(2*a); if (v1>0) newV = p + d * v1; else if (v2>0) newV = p + d * v2; else MITK_INFO << "ERROR1 - linear resampling"; j--; } //#pragma omp critical { vtkIdType id = vtkNewPoints->InsertNextPoint(newV.data_block()); container->GetPointIds()->InsertNextId(id); } lastV = newV; } else if ( (j==vertices.size()-1 && new_dist>0.0001) || seg_len<=0.0000001) { //#pragma omp critical { vtkIdType id = vtkNewPoints->InsertNextPoint(vertices.at(j).data_block()); container->GetPointIds()->InsertNextId(id); } } } //#pragma omp critical { newFiberWeights->SetValue(vtkNewCells->GetNumberOfCells(), weight); vtkNewCells->InsertNextCell(container); if (container->GetNumberOfPoints()!=targetPoints) unequal_fibs = true; } } // }); // MITK_INFO << "Done"; if (vtkNewCells->GetNumberOfCells()>0) { SetFiberWeights(newFiberWeights); m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } } } void mitk::FiberBundle::ResampleLinear(double pointDistance) { vtkSmartPointer vtkNewPoints = vtkSmartPointer::New(); vtkSmartPointer vtkNewCells = vtkSmartPointer::New(); MITK_INFO << "Resampling fibers (linear)"; boost::progress_display disp(static_cast(m_FiberPolyData->GetNumberOfCells())); vtkSmartPointer newFiberWeights = vtkSmartPointer::New(); newFiberWeights->SetName("FIBER_WEIGHTS"); newFiberWeights->SetNumberOfValues(m_NumFibers); std::vector< vtkSmartPointer > resampled_streamlines; resampled_streamlines.resize(static_cast(m_FiberPolyData->GetNumberOfCells())); #pragma omp parallel for for (int i=0; i(m_FiberPolyData->GetNumberOfCells()); i++) { std::vector< vnl_vector_fixed< double, 3 > > vertices; #pragma omp critical { ++disp; vtkCell* cell = m_FiberPolyData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); for (int j=0; jGetPoint(j, cand); vnl_vector_fixed< double, 3 > candV; candV[0]=cand[0]; candV[1]=cand[1]; candV[2]=cand[2]; vertices.push_back(candV); } } vtkSmartPointer container = vtkSmartPointer::New(); vnl_vector_fixed< double, 3 > lastV = vertices.at(0); #pragma omp critical { vtkIdType id = vtkNewPoints->InsertNextPoint(lastV.data_block()); container->GetPointIds()->InsertNextId(id); } for (unsigned int j=1; j vec = vertices.at(j) - lastV; double new_dist = vec.magnitude(); if (new_dist >= pointDistance) { vnl_vector_fixed< double, 3 > newV = lastV; if ( new_dist-pointDistance <= mitk::eps ) { vec.normalize(); newV += vec * pointDistance; } else { // intersection between sphere (radius 'pointDistance', center 'lastV') and line (direction 'd' and point 'p') vnl_vector_fixed< double, 3 > p = vertices.at(j-1); vnl_vector_fixed< double, 3 > d = vertices.at(j) - p; double a = d[0]*d[0] + d[1]*d[1] + d[2]*d[2]; double b = 2 * (d[0] * (p[0] - lastV[0]) + d[1] * (p[1] - lastV[1]) + d[2] * (p[2] - lastV[2])); double c = (p[0] - lastV[0])*(p[0] - lastV[0]) + (p[1] - lastV[1])*(p[1] - lastV[1]) + (p[2] - lastV[2])*(p[2] - lastV[2]) - pointDistance*pointDistance; double v1 =(-b + std::sqrt(b*b-4*a*c))/(2*a); double v2 =(-b - std::sqrt(b*b-4*a*c))/(2*a); if (v1>0) newV = p + d * v1; else if (v2>0) newV = p + d * v2; else MITK_INFO << "ERROR1 - linear resampling"; j--; } #pragma omp critical { vtkIdType id = vtkNewPoints->InsertNextPoint(newV.data_block()); container->GetPointIds()->InsertNextId(id); } lastV = newV; } else if (j==vertices.size()-1 && new_dist>0.0001) { #pragma omp critical { vtkIdType id = vtkNewPoints->InsertNextPoint(vertices.at(j).data_block()); container->GetPointIds()->InsertNextId(id); } } } #pragma omp critical { resampled_streamlines[static_cast(i)] = container; } } for (auto container : resampled_streamlines) { vtkNewCells->InsertNextCell(container); } if (vtkNewCells->GetNumberOfCells()>0) { m_FiberPolyData = vtkSmartPointer::New(); m_FiberPolyData->SetPoints(vtkNewPoints); m_FiberPolyData->SetLines(vtkNewCells); this->SetFiberPolyData(m_FiberPolyData, true); } } // reapply selected colorcoding in case PolyData structure has changed bool mitk::FiberBundle::Equals(mitk::FiberBundle* fib, double eps) { if (fib==nullptr) { MITK_INFO << "Reference bundle is nullptr!"; return false; } if (m_NumFibers!=fib->GetNumFibers()) { MITK_INFO << "Unequal number of fibers!"; MITK_INFO << m_NumFibers << " vs. " << fib->GetNumFibers(); return false; } for (unsigned int i=0; iGetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkCell* cell2 = fib->GetFiberPolyData()->GetCell(i); auto numPoints2 = cell2->GetNumberOfPoints(); vtkPoints* points2 = cell2->GetPoints(); if (numPoints2!=numPoints) { MITK_INFO << "Unequal number of points in fiber " << i << "!"; MITK_INFO << numPoints2 << " vs. " << numPoints; return false; } for (int j=0; jGetPoint(j); double* p2 = points2->GetPoint(j); if (fabs(p1[0]-p2[0])>eps || fabs(p1[1]-p2[1])>eps || fabs(p1[2]-p2[2])>eps) { MITK_INFO << "Unequal points in fiber " << i << " at position " << j << "!"; MITK_INFO << "p1: " << p1[0] << ", " << p1[1] << ", " << p1[2]; MITK_INFO << "p2: " << p2[0] << ", " << p2[1] << ", " << p2[2]; return false; } } } return true; } void mitk::FiberBundle::PrintSelf(std::ostream &os, itk::Indent indent) const { os << indent << "Number of fibers: " << this->GetNumFibers() << std::endl; os << indent << "Min. fiber length: " << this->GetMinFiberLength() << std::endl; os << indent << "Max. fiber length: " << this->GetMaxFiberLength() << std::endl; os << indent << "Mean fiber length: " << this->GetMeanFiberLength() << std::endl; os << indent << "Median fiber length: " << this->GetMedianFiberLength() << std::endl; os << indent << "STDEV fiber length: " << this->GetLengthStDev() << std::endl; os << indent << "Number of points: " << this->GetNumberOfPoints() << std::endl; os << indent << "Extent x: " << this->GetGeometry()->GetExtentInMM(0) << "mm" << std::endl; os << indent << "Extent y: " << this->GetGeometry()->GetExtentInMM(1) << "mm" << std::endl; os << indent << "Extent z: " << this->GetGeometry()->GetExtentInMM(2) << "mm" << std::endl; os << indent << "Diagonal: " << this->GetGeometry()->GetDiagonalLength() << "mm" << std::endl; os << "\nReference geometry:" << std::endl; os << indent << "Size: [" << std::defaultfloat << m_TrackVisHeader.dim[0] << " " << m_TrackVisHeader.dim[1] << " " << m_TrackVisHeader.dim[2] << "]" << std::endl; os << indent << "Voxel size: [" << m_TrackVisHeader.voxel_size[0] << " " << m_TrackVisHeader.voxel_size[1] << " " << m_TrackVisHeader.voxel_size[2] << "]" << std::endl; os << indent << "Origin: [" << m_TrackVisHeader.origin[0] << " " << m_TrackVisHeader.origin[1] << " " << m_TrackVisHeader.origin[2] << "]" << std::endl; os << indent << "Matrix: " << std::scientific << std::endl; os << indent << "[[" << m_TrackVisHeader.vox_to_ras[0][0] << ", " << m_TrackVisHeader.vox_to_ras[0][1] << ", " << m_TrackVisHeader.vox_to_ras[0][2] << ", " << m_TrackVisHeader.vox_to_ras[0][3] << "]" << std::endl; os << indent << " [" << m_TrackVisHeader.vox_to_ras[1][0] << ", " << m_TrackVisHeader.vox_to_ras[1][1] << ", " << m_TrackVisHeader.vox_to_ras[1][2] << ", " << m_TrackVisHeader.vox_to_ras[1][3] << "]" << std::endl; os << indent << " [" << m_TrackVisHeader.vox_to_ras[2][0] << ", " << m_TrackVisHeader.vox_to_ras[2][1] << ", " << m_TrackVisHeader.vox_to_ras[2][2] << ", " << m_TrackVisHeader.vox_to_ras[2][3] << "]" << std::endl; os << indent << " [" << m_TrackVisHeader.vox_to_ras[3][0] << ", " << m_TrackVisHeader.vox_to_ras[3][1] << ", " << m_TrackVisHeader.vox_to_ras[3][2] << ", " << m_TrackVisHeader.vox_to_ras[3][3] << "]]" << std::defaultfloat << std::endl; if (m_FiberWeights!=nullptr) { std::vector< float > weights; for (int i=0; iGetSize(); i++) weights.push_back(m_FiberWeights->GetValue(i)); std::sort(weights.begin(), weights.end()); os << "\nFiber weight statistics" << std::endl; os << indent << "Min: " << weights.front() << std::endl; os << indent << "1% quantile: " << weights.at(static_cast(weights.size()*0.01)) << std::endl; os << indent << "5% quantile: " << weights.at(static_cast(weights.size()*0.05)) << std::endl; os << indent << "25% quantile: " << weights.at(static_cast(weights.size()*0.25)) << std::endl; os << indent << "Median: " << weights.at(static_cast(weights.size()*0.5)) << std::endl; os << indent << "75% quantile: " << weights.at(static_cast(weights.size()*0.75)) << std::endl; os << indent << "95% quantile: " << weights.at(static_cast(weights.size()*0.95)) << std::endl; os << indent << "99% quantile: " << weights.at(static_cast(weights.size()*0.99)) << std::endl; os << indent << "Max: " << weights.back() << std::endl; } else os << indent << "\n\nNo fiber weight array found." << std::endl; Superclass::PrintSelf(os, 0); } mitk::FiberBundle::TrackVis_header mitk::FiberBundle::GetTrackVisHeader() { if (m_TrackVisHeader.hdr_size==0) { mitk::Geometry3D::Pointer geom = dynamic_cast(this->GetGeometry()); SetTrackVisHeader(geom); } return m_TrackVisHeader; } void mitk::FiberBundle::SetTrackVisHeader(const mitk::FiberBundle::TrackVis_header &TrackVisHeader) { m_TrackVisHeader = TrackVisHeader; } void mitk::FiberBundle::SetTrackVisHeader(mitk::BaseGeometry* geometry) { vtkSmartPointer< vtkMatrix4x4 > matrix = vtkSmartPointer< vtkMatrix4x4 >::New(); matrix->Identity(); if (geometry==nullptr) return; for(int i=0; i<3 ;i++) { m_TrackVisHeader.dim[i] = geometry->GetExtent(i); m_TrackVisHeader.voxel_size[i] = geometry->GetSpacing()[i]; m_TrackVisHeader.origin[i] = geometry->GetOrigin()[i]; matrix = geometry->GetVtkMatrix(); } for (int i=0; i<4; ++i) for (int j=0; j<4; ++j) m_TrackVisHeader.vox_to_ras[i][j] = matrix->GetElement(i, j); m_TrackVisHeader.n_scalars = 0; m_TrackVisHeader.n_properties = 0; sprintf(m_TrackVisHeader.voxel_order,"LPS"); m_TrackVisHeader.image_orientation_patient[0] = 1.0; m_TrackVisHeader.image_orientation_patient[1] = 0.0; m_TrackVisHeader.image_orientation_patient[2] = 0.0; m_TrackVisHeader.image_orientation_patient[3] = 0.0; m_TrackVisHeader.image_orientation_patient[4] = 1.0; m_TrackVisHeader.image_orientation_patient[5] = 0.0; m_TrackVisHeader.pad1[0] = 0; m_TrackVisHeader.pad1[1] = 0; m_TrackVisHeader.pad2[0] = 0; m_TrackVisHeader.pad2[1] = 0; m_TrackVisHeader.invert_x = 0; m_TrackVisHeader.invert_y = 0; m_TrackVisHeader.invert_z = 0; m_TrackVisHeader.swap_xy = 0; m_TrackVisHeader.swap_yz = 0; m_TrackVisHeader.swap_zx = 0; m_TrackVisHeader.n_count = 0; m_TrackVisHeader.version = 2; m_TrackVisHeader.hdr_size = 1000; std::string id = "TRACK"; strcpy(m_TrackVisHeader.id_string, id.c_str()); } /* ESSENTIAL IMPLEMENTATION OF SUPERCLASS METHODS */ void mitk::FiberBundle::UpdateOutputInformation() { } void mitk::FiberBundle::SetRequestedRegionToLargestPossibleRegion() { } bool mitk::FiberBundle::RequestedRegionIsOutsideOfTheBufferedRegion() { return false; } bool mitk::FiberBundle::VerifyRequestedRegion() { return true; } void mitk::FiberBundle::SetRequestedRegion(const itk::DataObject* ) { } diff --git a/Modules/FiberDissection/Interactor/mitkStreamlineInteractor.cpp b/Modules/FiberDissection/Interactor/mitkStreamlineInteractor.cpp index 2cfd0c4..ee82c54 100644 --- a/Modules/FiberDissection/Interactor/mitkStreamlineInteractor.cpp +++ b/Modules/FiberDissection/Interactor/mitkStreamlineInteractor.cpp @@ -1,1185 +1,1205 @@ /*============================================================================ The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center (DKFZ) All rights reserved. Use of this source code is governed by a 3-clause BSD license that can be found in the LICENSE file. ============================================================================*/ #include "mitkStreamlineInteractor.h" //#include "mitkStreamlineMapper2D.h" // MITK includes #include #include #include #include #include #include #include #include #include #include // VTK includes #include #include #include #include #include #include #include #include #include #include mitk::StreamlineInteractor::StreamlineInteractor() { m_ColorForHighlight[0] = 1.0; m_ColorForHighlight[1] = 0.5; m_ColorForHighlight[2] = 0.0; m_ColorForHighlight[3] = 1.0; // TODO if we want to get this configurable, the this is the recipe: // - make the 2D mapper add corresponding properties to control "enabled" and "color" // - make the interactor evaluate those properties // - in an ideal world, modify the state machine on the fly and skip mouse move handling } mitk::StreamlineInteractor::~StreamlineInteractor() { } void mitk::StreamlineInteractor::ConnectActionsAndFunctions() { // CONNECT_CONDITION("isoverstreamline", HasPickedHandle); CONNECT_CONDITION("isoverstreamline", CheckSelection); CONNECT_FUNCTION("selectstreamline", SelectStreamline); CONNECT_FUNCTION("addnegstreamline", AddStreamlineNegBundle); CONNECT_FUNCTION("addposstreamline", AddStreamlinePosBundle); CONNECT_FUNCTION("addnegtolabelstreamline", AddNegStreamlinetolabelsBundle); CONNECT_FUNCTION("addpostolabelstreamline", AddPosStreamlinetolabelsBundle); // CONNECT_FUNCTION("FeedUndoStack", FeedUndoStack); } void mitk::StreamlineInteractor::SetNegativeNode(DataNode *node) { // DataInteractor::SetDataNode(node); m_NegStreamlineNode = node; m_NegStreamline= dynamic_cast(node->GetData()); MITK_INFO << "Negative Node added"; } void mitk::StreamlineInteractor::SetToLabelNode(DataNode *node) { m_manStreamlineNode = node; DataInteractor::SetDataNode(m_manStreamlineNode); m_manStreamline = dynamic_cast(node->GetData()); MITK_INFO << "Label node added"; +// m_init = false; } void mitk::StreamlineInteractor::SetPositiveNode(DataNode *node) { // DataInteractor::SetDataNode(node); m_PosStreamlineNode = node; m_PosStreamline= dynamic_cast(node->GetData()); MITK_INFO << "Positive Node added"; - } +} + +void mitk::StreamlineInteractor::LabelfromPrediction(bool predlabeling) +{ + m_predlabeling = predlabeling; + +} bool mitk::StreamlineInteractor::CheckSelection(const InteractionEvent *interactionEvent) { /* Save Colorvector here*/ - vtkSmartPointer color = m_manStreamline->GetFiberColors(); const auto *positionEvent = dynamic_cast(interactionEvent); if (positionEvent != nullptr) { if (interactionEvent->GetSender()->GetMapperID() == BaseRenderer::Standard2D) { BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard2D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); - if (picker->GetCellId()==-1) + if (picker->GetCellId()==-1 && m_predlabeling==false) { - m_manStreamline->SetFiberColors(255, 255, 255); +// m_manStreamline->SetFiberColors(255, 255, 255); +// vtkSmartPointer FiberColors = m_manStreamline->m_FiberColors; + RenderingManager::GetInstance()->RequestUpdateAll(); return false; } else { return true; } } else { BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard3D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); - if (picker->GetCellId()==-1) + if (picker->GetCellId()==-1 && m_predlabeling==false) { - m_manStreamline->SetFiberColors(255, 255, 255); +// m_manStreamline->SetFiberColors(255, 255, 255); RenderingManager::GetInstance()->RequestUpdateAll(); return false; } else { return true; } } } else { // m_manStreamline->SetFiberColors(255, 255, 255); RenderingManager::GetInstance()->RequestUpdateAll(); return false; } } void mitk::StreamlineInteractor::SelectStreamline(StateMachineAction *, InteractionEvent *interactionEvent) { auto *positionEvent = dynamic_cast(interactionEvent); if (interactionEvent->GetSender()->GetMapperID() == BaseRenderer::Standard2D) { BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard2D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); - + if (m_predlabeling==false) + { vtkIdType pickedCellID = picker->GetCellId(); - m_manStreamline->SetSingleFiberColor(0, 255, 0, pickedCellID); + m_manStreamline->SetSingleFiberColor(0.0, 255.0, 0.0, pickedCellID); + } } else { BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard3D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); - + if (m_predlabeling==false) + { vtkIdType pickedCellID = picker->GetCellId(); - - m_manStreamline->SetSingleFiberColor(0, 255, 0, pickedCellID); + m_manStreamline->SetSingleFiberColor(0.0, 255.0, 0.0, pickedCellID); + } } + RenderingManager::GetInstance()->RequestUpdateAll(); } void mitk::StreamlineInteractor::AddStreamlinePosBundle(StateMachineAction *, InteractionEvent *interactionEvent) { MITK_INFO << "PositiveBundle clicked"; DataInteractor::SetDataNode(m_manStreamlineNode); auto positionEvent = dynamic_cast(interactionEvent); if (positionEvent == nullptr) { MITK_INFO << "no position"; } if (interactionEvent->GetSender()->GetMapperID() == BaseRenderer::Standard2D) { // m_PickedHandle = PickFrom2D(positionEvent); MITK_INFO << "2D"; BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard2D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); vtkIdType pickedCellID = picker->GetCellId(); if (picker->GetCellId()==-1) { MITK_INFO << "Nothing picked"; } else { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_PosStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(pickedCellID); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // m_PosStreamline = mitk::FiberBundle::New(vNewPolyData); m_PosStreamline->GetFiberPolyData()->SetPoints(vNewPoints); m_PosStreamline->GetFiberPolyData()->SetLines(vNewLines); m_PosStreamline->SetFiberColors(0, 255, 0); // m_PosStreamline->SetFiberWeights(m_PosStreamline->GetFiberWeights()); m_manStreamline->GetFiberPolyData()->DeleteCell(pickedCellID); m_manStreamline->GetFiberPolyData()->RemoveDeletedCells(); - + if (m_predlabeling==false) + { vtkSmartPointer vNewPolyData2 = vtkSmartPointer::New(); vtkSmartPointer vNewLines2 = vtkSmartPointer::New(); vtkSmartPointer vNewPoints2 = vtkSmartPointer::New(); counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints2->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines2->InsertNextCell(container); counter++; } m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints2); m_manStreamline->GetFiberPolyData()->SetLines(vNewLines2); m_manStreamline->SetFiberColors(255, 255, 255); + } } } else { BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard3D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); vtkIdType pickedCellID = picker->GetCellId(); if (picker->GetCellId()==-1) { MITK_INFO << "Nothing picked"; } else { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_PosStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(pickedCellID); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // m_PosStreamline = mitk::FiberBundle::New(vNewPolyData); m_PosStreamline->GetFiberPolyData()->SetPoints(vNewPoints); m_PosStreamline->GetFiberPolyData()->SetLines(vNewLines); m_PosStreamline->SetFiberColors(0, 255, 0); // m_PosStreamline->SetFiberWeights(m_PosStreamline->GetFiberWeights()); m_manStreamline->GetFiberPolyData()->DeleteCell(pickedCellID); m_manStreamline->GetFiberPolyData()->RemoveDeletedCells(); - + if (m_predlabeling==false) + { vtkSmartPointer vNewPolyData2 = vtkSmartPointer::New(); vtkSmartPointer vNewLines2 = vtkSmartPointer::New(); vtkSmartPointer vNewPoints2 = vtkSmartPointer::New(); counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints2->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines2->InsertNextCell(container); counter++; } m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints2); m_manStreamline->GetFiberPolyData()->SetLines(vNewLines2); m_manStreamline->SetFiberColors(255, 255, 255); + } } } RenderingManager::GetInstance()->RequestUpdateAll(); } void mitk::StreamlineInteractor::AddStreamlineNegBundle(StateMachineAction *, InteractionEvent *interactionEvent) { MITK_INFO << "NegativeBundle clicked"; DataInteractor::SetDataNode(m_manStreamlineNode); auto positionEvent = dynamic_cast(interactionEvent); if (positionEvent == nullptr) { MITK_INFO << "no position"; } if (interactionEvent->GetSender()->GetMapperID() == BaseRenderer::Standard2D) { // m_PickedHandle = PickFrom2D(positionEvent); MITK_INFO << "2D"; BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard2D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); vtkIdType pickedCellID = picker->GetCellId(); if (picker->GetCellId()==-1) { MITK_INFO << "Nothing picked"; } else { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_NegStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(pickedCellID); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // m_NegStreamline = mitk::FiberBundle::New(vNewPolyData); m_NegStreamline->GetFiberPolyData()->SetPoints(vNewPoints); m_NegStreamline->GetFiberPolyData()->SetLines(vNewLines); m_NegStreamline->SetFiberColors(255, 0, 0); // m_NegStreamline->SetFiberWeights(m_NegStreamline->GetFiberWeights()); m_manStreamline->GetFiberPolyData()->DeleteCell(pickedCellID); m_manStreamline->GetFiberPolyData()->RemoveDeletedCells(); + if (m_predlabeling==false) + { vtkSmartPointer vNewPolyData2 = vtkSmartPointer::New(); vtkSmartPointer vNewLines2 = vtkSmartPointer::New(); vtkSmartPointer vNewPoints2 = vtkSmartPointer::New(); counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints2->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines2->InsertNextCell(container); counter++; } m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints2); m_manStreamline->GetFiberPolyData()->SetLines(vNewLines2); m_manStreamline->SetFiberColors(255, 255, 255); - + } } } else { BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; // if (picker == nullptr) // { picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard3D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); vtkIdType pickedCellID = picker->GetCellId(); if (picker->GetCellId()==-1) { MITK_INFO << "Nothing picked"; } else { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_NegStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(pickedCellID); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // m_NegStreamline = mitk::FiberBundle::New(vNewPolyData); m_NegStreamline->GetFiberPolyData()->SetPoints(vNewPoints); m_NegStreamline->GetFiberPolyData()->SetLines(vNewLines); m_NegStreamline->SetFiberColors(255, 0, 0); m_manStreamline->GetFiberPolyData()->DeleteCell(pickedCellID); m_manStreamline->GetFiberPolyData()->RemoveDeletedCells(); - vtkSmartPointer vNewPolyData2 = vtkSmartPointer::New(); - vtkSmartPointer vNewLines2 = vtkSmartPointer::New(); - vtkSmartPointer vNewPoints2 = vtkSmartPointer::New(); - - - counter = 0; - for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) - { - - vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); - auto numPoints = cell->GetNumberOfPoints(); - vtkPoints* points = cell->GetPoints(); - - vtkSmartPointer container = vtkSmartPointer::New(); - for (unsigned int j=0; jGetPoint(j, p); - - vtkIdType id = vNewPoints2->InsertNextPoint(p); - container->GetPointIds()->InsertNextId(id); + vtkSmartPointer vNewPolyData2 = vtkSmartPointer::New(); + vtkSmartPointer vNewLines2 = vtkSmartPointer::New(); + vtkSmartPointer vNewPoints2 = vtkSmartPointer::New(); + + + counter = 0; + for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) + { + + vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); + auto numPoints = cell->GetNumberOfPoints(); + vtkPoints* points = cell->GetPoints(); + + vtkSmartPointer container = vtkSmartPointer::New(); + for (unsigned int j=0; jGetPoint(j, p); + + vtkIdType id = vNewPoints2->InsertNextPoint(p); + container->GetPointIds()->InsertNextId(id); + } + // weights->InsertValue(counter, fib->GetFiberWeight(i)); + vNewLines2->InsertNextCell(container); + counter++; + + } + m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints2); + m_manStreamline->GetFiberPolyData()->SetLines(vNewLines2); + m_manStreamline->SetFiberColors(255, 255, 255); } - // weights->InsertValue(counter, fib->GetFiberWeight(i)); - vNewLines2->InsertNextCell(container); - counter++; - - } - m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints2); - m_manStreamline->GetFiberPolyData()->SetLines(vNewLines2); - m_manStreamline->SetFiberColors(255, 255, 255); - } } RenderingManager::GetInstance()->RequestUpdateAll(); } void mitk::StreamlineInteractor::AddNegStreamlinetolabelsBundle(StateMachineAction *, InteractionEvent *interactionEvent) { MITK_INFO << "TolabelBundle clicked"; DataInteractor::SetDataNode(m_NegStreamlineNode); auto positionEvent = dynamic_cast(interactionEvent); if (positionEvent == nullptr) { MITK_INFO << "no position"; } if (interactionEvent->GetSender()->GetMapperID() == BaseRenderer::Standard2D) { // m_PickedHandle = PickFrom2D(positionEvent); MITK_INFO << "2D"; BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard2D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); vtkIdType pickedCellID = picker->GetCellId(); if (picker->GetCellId()==-1) { MITK_INFO << "Nothing picked"; } else { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vtkCell* cell = m_NegStreamline->GetFiberPolyData()->GetCell(pickedCellID); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // m_manStreamline = mitk::FiberBundle::New(vNewPolyData); m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints); m_manStreamline->GetFiberPolyData()->SetLines(vNewLines); m_manStreamline->SetFiberColors(255, 255, 255); // m_manStreamline->SetFiberWeights(m_manStreamline->GetFiberWeights()); m_NegStreamline->GetFiberPolyData()->DeleteCell(pickedCellID); m_NegStreamline->GetFiberPolyData()->RemoveDeletedCells(); } } else { BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; // if (picker == nullptr) // { picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard3D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); vtkIdType pickedCellID = picker->GetCellId(); if (picker->GetCellId()==-1) { MITK_INFO << "Nothing picked"; } else { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vtkCell* cell = m_NegStreamline->GetFiberPolyData()->GetCell(pickedCellID); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // m_manStreamline = mitk::FiberBundle::New(vNewPolyData); m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints); m_manStreamline->GetFiberPolyData()->SetLines(vNewLines); m_manStreamline->SetFiberColors(255, 255, 255); m_NegStreamline->GetFiberPolyData()->DeleteCell(pickedCellID); m_NegStreamline->GetFiberPolyData()->RemoveDeletedCells(); } } RenderingManager::GetInstance()->RequestUpdateAll(); } void mitk::StreamlineInteractor::AddPosStreamlinetolabelsBundle(StateMachineAction *, InteractionEvent *interactionEvent) { MITK_INFO << "TolabelBundle clicked"; DataInteractor::SetDataNode(m_PosStreamlineNode); auto positionEvent = dynamic_cast(interactionEvent); if (positionEvent == nullptr) { MITK_INFO << "no position"; } if (interactionEvent->GetSender()->GetMapperID() == BaseRenderer::Standard2D) { // m_PickedHandle = PickFrom2D(positionEvent); MITK_INFO << "2D"; BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard2D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); vtkIdType pickedCellID = picker->GetCellId(); if (picker->GetCellId()==-1) { MITK_INFO << "Nothing picked"; } else { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vtkCell* cell = m_PosStreamline->GetFiberPolyData()->GetCell(pickedCellID); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // m_manStreamline = mitk::FiberBundle::New(vNewPolyData); m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints); m_manStreamline->GetFiberPolyData()->SetLines(vNewLines); m_manStreamline->SetFiberColors(255, 255, 255); // m_manStreamline->SetFiberWeights(m_manStreamline->GetFiberWeights()); m_PosStreamline->GetFiberPolyData()->DeleteCell(pickedCellID); m_PosStreamline->GetFiberPolyData()->RemoveDeletedCells(); } } else { BaseRenderer *renderer = positionEvent->GetSender(); auto &picker = m_Picker[renderer]; // if (picker == nullptr) // { picker = vtkSmartPointer::New(); picker->SetTolerance(0.01); auto mapper = GetDataNode()->GetMapper(BaseRenderer::Standard3D); auto vtk_mapper = dynamic_cast(mapper); if (vtk_mapper) { // doing this each time is bizarre picker->AddPickList(vtk_mapper->GetVtkProp(renderer)); picker->PickFromListOn(); } // } auto displayPosition = positionEvent->GetPointerPositionOnScreen(); picker->Pick(displayPosition[0], displayPosition[1], 0, positionEvent->GetSender()->GetVtkRenderer()); vtkIdType pickedCellID = picker->GetCellId(); if (picker->GetCellId()==-1) { MITK_INFO << "Nothing picked"; } else { vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for ( int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = m_manStreamline->GetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vtkCell* cell = m_PosStreamline->GetFiberPolyData()->GetCell(pickedCellID); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); // m_manStreamline = mitk::FiberBundle::New(vNewPolyData); m_manStreamline->GetFiberPolyData()->SetPoints(vNewPoints); m_manStreamline->GetFiberPolyData()->SetLines(vNewLines); m_manStreamline->SetFiberColors(255, 255, 255); m_PosStreamline->GetFiberPolyData()->DeleteCell(pickedCellID); m_PosStreamline->GetFiberPolyData()->RemoveDeletedCells(); } } RenderingManager::GetInstance()->RequestUpdateAll(); } diff --git a/Modules/FiberDissection/Interactor/mitkStreamlineInteractor.h b/Modules/FiberDissection/Interactor/mitkStreamlineInteractor.h index 23f0595..3b43f1e 100644 --- a/Modules/FiberDissection/Interactor/mitkStreamlineInteractor.h +++ b/Modules/FiberDissection/Interactor/mitkStreamlineInteractor.h @@ -1,118 +1,123 @@ /*============================================================================ The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center (DKFZ) All rights reserved. Use of this source code is governed by a 3-clause BSD license that can be found in the LICENSE file. ============================================================================*/ #ifndef mitkStreamlineInteractor_h #define mitkStreamlineInteractor_h #include "MitkFiberDissectionExports.h" // MITK includes #include #include #include #include // VTK includes #include #include // System includes #include namespace mitk { class InteractionPositionEvent; //! Data interactor to pick streamlines via interaction //! with a mitk::Streamline. //! //! //! To determine what parts of the object are clicked during interaction, //! the mappers (2D: custom mapper, 3D: regular surface mapper) are asked //! for their VTK objects, picking is performed, and the picked point is //! forwarded to the Streamline object for interpretation. //! //! The interactor fills the undo/redo stack with operations on the modified geometry. //! //! \sa Streamline class MITKFIBERDISSECTION_EXPORT StreamlineInteractor : public DataInteractor { public: mitkClassMacro(StreamlineInteractor, DataInteractor); itkFactorylessNewMacro(Self); itkCloneMacro(Self); //! The node holding the Fiberbundle for visual feedback. //! This is the node that the interactor is primarily working on //! (calls DataInteractor::SetDataNode). void SetNegativeNode(DataNode *node); void SetToLabelNode(DataNode *node); void SetPositiveNode(DataNode *node); + void LabelfromPrediction(bool predlabeling); protected: void AddStreamlineNegBundle(StateMachineAction *, InteractionEvent *interactionEvent); void AddStreamlinePosBundle(StateMachineAction *, InteractionEvent *interactionEvent); void AddNegStreamlinetolabelsBundle(StateMachineAction *, InteractionEvent *interactionEvent); void AddPosStreamlinetolabelsBundle(StateMachineAction *, InteractionEvent *interactionEvent); bool CheckSelection(const InteractionEvent *interactionEvent); void SelectStreamline(StateMachineAction *, InteractionEvent *); + std::map> m_Picker; private: StreamlineInteractor(); ~StreamlineInteractor() override; //! Setup the relation between the XML state machine and this object's methods. void ConnectActionsAndFunctions() override; //! State machine condition: successful Streamline picking //! \return true when any part of the Streamline has been picked. // bool HasPickedHandle(const InteractionEvent *); // void DecideInteraction(StateMachineAction *, InteractionEvent *interactionEvent); // //! Pick a Streamline handle from a 2D event (passing by the 2D mapper) // Streamline::HandleType PickFrom2D(const InteractionPositionEvent *positionEvent); // //! Pick a Streamline handle from a 3D event // //! (passing by the general surface mapper and the Streamline object) // Streamline::HandleType PickFrom3D(const InteractionPositionEvent *positionEvent); // void UpdateHandleHighlight(); //! the Streamline used for visual feedback and picking mitk::FiberBundle::Pointer m_NegStreamline; mitk::FiberBundle::Pointer m_PosStreamline; mitk::FiberBundle::Pointer m_manStreamline; + bool m_predlabeling; // mitk::FiberBundle::Pointer m_visualStreamline; mitk::DataNode::Pointer m_NegStreamlineNode; mitk::DataNode::Pointer m_PosStreamlineNode; mitk::DataNode::Pointer m_manStreamlineNode; // mitk::DataNode::Pointer m_visualStreamlineNode; + + vtkSmartPointer m_extracted_streamline; double m_ColorForHighlight[4]; }; } #endif diff --git a/Modules/FiberDissection/MachineLearning/mitkStreamlineFeatureExtractor.cpp b/Modules/FiberDissection/MachineLearning/mitkStreamlineFeatureExtractor.cpp index b02ebed..10707e4 100644 --- a/Modules/FiberDissection/MachineLearning/mitkStreamlineFeatureExtractor.cpp +++ b/Modules/FiberDissection/MachineLearning/mitkStreamlineFeatureExtractor.cpp @@ -1,767 +1,799 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #include "mitkStreamlineFeatureExtractor.h" #define _USE_MATH_DEFINES #include #include #include #include namespace mitk{ StreamlineFeatureExtractor::StreamlineFeatureExtractor() : m_NumPoints(40) { } StreamlineFeatureExtractor::~StreamlineFeatureExtractor() { } void StreamlineFeatureExtractor::SetTractogramPlus(const mitk::FiberBundle::Pointer &TractogramPlus) { m_TractogramPlus = TractogramPlus; } void StreamlineFeatureExtractor::SetTractogramMinus(const mitk::FiberBundle::Pointer &TractogramMinus) { m_TractogramMinus = TractogramMinus; } void StreamlineFeatureExtractor::SetTractogramPrototypes(const mitk::FiberBundle::Pointer &TractogramPrototypes, bool standard) { if (standard) { MITK_INFO << "Use Standard Prototypes..."; - m_inputPrototypes = mitk::IOUtil::Load("/home/r948e/E132-Projekte/Projects/2022_Peretzke_Interactive_Fiber_Dissection/mitk_diff/prototypes_599671.trk"); + m_inputPrototypes = mitk::IOUtil::Load("/home/r948e/E132-Projekte/Projects/2022_Peretzke_Interactive_Fiber_Dissection/mitk_diff/prototypes_599671_40.trk"); } else { MITK_INFO << "Use individual Prototypes..."; m_inputPrototypes = TractogramPrototypes; } } void StreamlineFeatureExtractor::SetActiveCycle(int &activeCycle) { m_activeCycle= activeCycle; } +void StreamlineFeatureExtractor::SetInitRandom(int &initRandom) +{ + m_initRandom= initRandom; +} + void StreamlineFeatureExtractor::SetTractogramTest(const mitk::FiberBundle::Pointer &TractogramTest, std::string TractogramTestName) { std::string path = "/home/r948e/E132-Projekte/Projects/2022_Peretzke_Interactive_Fiber_Dissection/mitk_diff/storage/"; path.append(TractogramTestName); m_TractogramTest= TractogramTest; auto s = std::to_string(m_NumPoints); - m_DistancesTestName= path.append("_distances" + s + ".csv"); + m_DistancesTestName= path.append("_distances" + std::to_string(m_NumPoints) + "_" + std::to_string(m_activeCycle) + ".csv"); } std::vector > StreamlineFeatureExtractor::ResampleFibers(mitk::FiberBundle::Pointer tractogram) { mitk::FiberBundle::Pointer temp_fib = tractogram->GetDeepCopy(); temp_fib->ResampleToNumPoints(m_NumPoints); MITK_INFO << "Resampling Done"; std::vector< vnl_matrix > out_fib(temp_fib->GetFiberPolyData()->GetNumberOfCells()); for (int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = temp_fib->GetFiberPolyData()->GetCell(i); int numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vnl_matrix streamline; streamline.set_size(3, m_NumPoints); streamline.fill(0.0); for (int j=0; jGetPoint(j, cand); vnl_vector_fixed< float, 3 > candV; candV[0]=cand[0]; candV[1]=cand[1]; candV[2]=cand[2]; streamline.set_column(j, candV); } // out_fib.push_back(streamline); out_fib.at(i)=streamline; } // }); return out_fib; } std::vector > StreamlineFeatureExtractor::CalculateDmdf(std::vector > tractogram, std::vector > prototypes) { std::vector< vnl_matrix > dist_vec(tractogram.size());// MITK_INFO << "Start Calculating Dmdf"; cv::parallel_for_(cv::Range(0, tractogram.size()), [&](const cv::Range &range) { for (int i = range.start; i < range.end; i++) // for (unsigned int i=0; i distances; distances.set_size(1, prototypes.size()); distances.fill(0.0); for (unsigned int j=0; j single_distances; single_distances.set_size(1, tractogram.at(0).cols()); single_distances.fill(0.0); vnl_matrix single_distances_flip; single_distances_flip.set_size(1, tractogram.at(0).cols()); single_distances_flip.fill(0.0); for (unsigned int ik=0; ik single_distances.mean()) { distances.put(0,j, single_distances.mean()); } else { distances.put(0,j, single_distances_flip.mean()); } } // dist_vec.push_back(distances); dist_vec.at(i) = distances; } }); MITK_INFO << "Done Calculation"; MITK_INFO << dist_vec.at(0).size(); return dist_vec; } std::vector > StreamlineFeatureExtractor::MergeTractogram(std::vector > prototypes, std::vector > positive_local_prototypes, std::vector > negative_local_prototypes) { unsigned int pos_locals; unsigned int neg_locals; if (positive_local_prototypes.size() >= 50) { pos_locals= 50; } else { pos_locals= positive_local_prototypes.size(); } if (pos_locals <= negative_local_prototypes.size()) { neg_locals = pos_locals; } else { neg_locals= negative_local_prototypes.size(); } std::vector< vnl_matrix > merged_prototypes; for (unsigned int k=0; k StreamlineFeatureExtractor::Sort(std::vector sortingVector, int lengths) { std::vector index; std::priority_queue> q; for (unsigned int i = 0; i < sortingVector.size(); ++i) { q.push(std::pair(sortingVector[i], i)); } for (int i = 0; i < lengths; ++i) { int ki = q.top().second; index.push_back(ki); q.pop(); } return index; } std::vector> StreamlineFeatureExtractor::GetData() { MITK_INFO << "Start Function Get Data"; /*Vector which saves Prediction and Fibers to label based on uncertainty*/ std::vector> index_vec; // int labels_arr [m_DistancesPlus.size()+m_DistancesMinus.size()]; cv::Mat data; cv::Mat labels_arr_vec; int size_plus = 0; /*Create Trainingdata: Go through positive and negative Bundle and save distances as cv::Mat and create vector with labels*/ for ( unsigned int i=0; i seeds; for (int cont = 0; cont < labels_arr_vec.rows; cont++) { seeds.push_back(cont); } cv::randShuffle(seeds); cv::Mat labels_shuffled; cv::Mat samples_shuffled; + + + for (int cont = 0; cont < labels_arr_vec.rows; cont++) { labels_shuffled.push_back(labels_arr_vec.row(seeds[cont])); } for (int cont = 0; cont < labels_arr_vec.rows; cont++) { samples_shuffled.push_back(data.row(seeds[cont])); } + std::ofstream labelsfile; + labelsfile.open("/home/r948e/E132-Projekte/Projects/2022_Peretzke_Interactive_Fiber_Dissection/mitk_diff/storage/Labels_" + std::to_string(m_activeCycle) + ".csv"); + labelsfile<< cv::format(labels_shuffled, cv::Formatter::FMT_CSV) << std::endl; + labelsfile.close(); + + std::ofstream featuresfile; + featuresfile.open("/home/r948e/E132-Projekte/Projects/2022_Peretzke_Interactive_Fiber_Dissection/mitk_diff/storage/Features_" + std::to_string(m_activeCycle) + ".csv"); + featuresfile<< cv::format(samples_shuffled, cv::Formatter::FMT_CSV) << std::endl; + featuresfile.close(); /*Create Dataset and initialize Classifier*/ cv::Ptr m_traindata = cv::ml::TrainData::create(samples_shuffled, cv::ml::ROW_SAMPLE, labels_shuffled); auto statistic_model = cv::ml::RTrees::create(); auto criteria = cv::TermCriteria(); criteria.type = cv::TermCriteria::MAX_ITER; // criteria.epsilon = 1e-8; criteria.maxCount = 800; statistic_model->setMaxDepth(50); //set to three // statistic_model->setMinSampleCount(m_traindata->getNTrainSamples()*0.01); statistic_model->setMinSampleCount(2); statistic_model->setTruncatePrunedTree(false); statistic_model->setUse1SERule(false); statistic_model->setUseSurrogates(false); statistic_model->setTermCriteria(criteria); statistic_model->setCVFolds(1); statistic_model->setPriors(newweight); /*Train Classifier*/ MITK_INFO << "Start Training"; statistic_model->train(m_traindata); /*Predict on Test Data*/ MITK_INFO << "Predicting"; /*Create Dataset as cv::Mat*/ cv::Mat dataTest; for ( unsigned int i=0; i indexPrediction; std::vector e(m_DistancesTest.size()); /*For every Sample/Streamline get Prediction and entropy (=based on counts of Random Forest)*/ MITK_INFO << "Predicting on all cores"; #pragma omp parallel for for (unsigned int i=0; ipredict(dataTest.row(i)); #pragma omp critical if (val==1) { indexPrediction.push_back(i); } cv::Mat vote; statistic_model->getVotes(dataTest.row(i), vote, 0); e.at(i) = ( -(vote.at(1,0)*1.0)/ (vote.at(1,0)+vote.at(1,1)) * log2((vote.at(1,0)*1.0)/ (vote.at(1,0)+vote.at(1,1))) - (vote.at(1,1)*1.0)/ (vote.at(1,0)+vote.at(1,1))* log2((vote.at(1,1)*1.0)/ (vote.at(1,0)+vote.at(1,1)))); if (isnan(e.at(i))) { e.at(i)=0; } } MITK_INFO << "Done"; /*Save entropy values for analysis*/ -// std::ofstream myfile3; -// myfile3.open("/home/r948e/mycsv/entropydata" + std::to_string(m_activeCycle) + ".csv"); -// for (unsigned int i = 0; i < e.size(); i++) -// { -// myfile3 << e.at(i) << ' '; -// } -// myfile3.close(); + std::ofstream myfile3; + myfile3.open("/home/r948e/mycsv/entropydata_" + std::to_string(m_activeCycle) + ".csv"); + for (unsigned int i = 0; i < e.size(); i++) + { + myfile3 << e.at(i) << ' '; + } + myfile3.close(); MITK_INFO << "--------------"; MITK_INFO << "Prediction vector size:"; MITK_INFO << indexPrediction.size(); MITK_INFO << "Entropy vector size:"; MITK_INFO << e.size(); MITK_INFO << "--------------"; /*Get index of most unertain data (lengths defines how many data is saved)*/ - int lengths=500; -// int lengths = std::count_if(e.begin(), e.end(),[&](auto const& val){ return val >= 0.9; }); +// int lengths=500; + int lengths = std::count_if(e.begin(), e.end(),[&](auto const& val){ return val >= 0.95; }); + if (length<500) + { + length=500 + } + std::vector indexUnc = Sort(e, lengths); MITK_INFO << indexUnc.size(); vnl_matrix distances_matrix; distances_matrix.set_size(lengths, lengths); distances_matrix.fill(0.0); std::vector distances_matrix_mean; for (int i=0; i diff = m_DistancesTest.at(indexUnc.at(i)) - m_DistancesTest.at(indexUnc.at(k)); /*Into the eucledean difference matrix, put the distance in Feature Space between every sample pair*/ distances_matrix.put(i,k,diff.absolute_value_sum()/m_DistancesTest.at(0).size()); } /*For every Sample/Streamline get the mean eucledean distance to all other Samples => one value for every Sample*/ distances_matrix_mean.push_back(distances_matrix.get_row(i).mean()); // MITK_INFO << meanval.at(i); } // vnl_vector sum_matrix; // sum_matrix.set_size(lengths); // sum_matrix.fill(0.0); // MITK_INFO << distances_matrix.get_row(myidx.at(i)+ sum_matrix.get_row(0) // for (int i=0; i sum_matrix; // sum_matrix.set_size(lengths, lengths); // sum_matrix.set_size(lengths, lengths); /*Index to find values in distancematrix*/ std::vector myidx; /*Index to find actual streamlines using indexUnc*/ std::vector indexUncDist; /*Start with the Streamline of the highest entropy, which is in distance_matrix at idx 0*/ myidx.push_back(0); indexUncDist.push_back(indexUnc.at(myidx.at(0))); /*Vecotr that stores minvalues of current iteration*/ vnl_matrix cur_vec; cur_vec.set_size(1,lengths); cur_vec.fill(0.0); for (int i=0; i sum_matrix; sum_matrix.set_size(myidx.size(), lengths); sum_matrix.fill(0); for (unsigned int ii=0; ii indexUncdist; // std::priority_queue> qq; // for (unsigned int i = 0; i < distances_matrix_mean.size(); ++i) // { // qq.push(std::pair(distances_matrix_mean[i], i)); // } //// int k = m_DistancesTest.size(); // number of indices we need //// int k = lengths; // number of indices we need // for (int i = 0; i < lengths; ++i) // { // int kki = qq.top().second; // indexUncdist.push_back(indexUnc.at(kki)); // qq.pop(); // } MITK_INFO << "Dist_stop"; // std::ofstream myfile6; // myfile6.open("/home/r948e/mycsv/distances_matrix_mean.csv"); // for (unsigned int i = 0; i < distances_matrix_mean.size(); i++) // { // myfile6 << distances_matrix_mean.at(i) << ' '; // } // myfile6.close(); // std::ofstream myfile5; // myfile5.open("/home/r948e/mycsv/indexUncdist.csv"); // for (unsigned int i = 0; i < indexUncdist.size(); i++) // { // myfile5 << indexUncdist.at(i) << ' '; // } // myfile5.close(); // MITK_INFO << distances_matrix; // MITK_INFO << distances_matrix.max_value(); // MITK_INFO << distances_matrix.arg_max(); // vnl_matrix myx = m_DistancesTest.at(indexUnc.at(0)) - m_DistancesTest.at(indexUnc.at(1)); //// myx = (m_DistancesTest.at(indexUnc.at(0)) - m_DistancesTest.at(indexUnc.at(1))); // MITK_INFO << m_DistancesTest.at(indexUnc.at(0)); // MITK_INFO << m_DistancesTest.at(indexUnc.at(1)); // MITK_INFO << myx.get(0,0); // MITK_INFO << sqrt(pow(myx.get(0,0),2)); // MITK_INFO << myx.get(0,1); // MITK_INFO << myx.absolute_value_sum()/m_DistancesTest.at(0).size(); // MITK_INFO << "Done"; /*Save Prediction*/ index_vec.push_back(indexPrediction); /*Save index of uncertainty measures*/ index_vec.push_back(indexUnc); /*Save index of uncertainty measures influenced by distance*/ index_vec.push_back(indexUncDist); return index_vec; } mitk::FiberBundle::Pointer StreamlineFeatureExtractor::CreatePrediction(std::vector &index) { mitk::FiberBundle::Pointer Prediction; MITK_INFO << "Create Bundle"; vtkSmartPointer FibersData; FibersData = vtkSmartPointer::New(); FibersData->SetPoints(vtkSmartPointer::New()); FibersData->SetLines(vtkSmartPointer::New()); vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); // vtkSmartPointer weights = vtkSmartPointer::New(); // weights->SetNumberOfValues(this->GetNumFibers()+fib->GetNumFibers()); unsigned int indexSize = index.size(); unsigned int counter = 0; MITK_INFO << "Start Loop"; for (unsigned int i=0; iGetFiberPolyData()->GetCell(index[i]); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } // weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } MITK_INFO << "Counter"; MITK_INFO << counter; vNewPolyData->SetLines(vNewLines); vNewPolyData->SetPoints(vNewPoints); FibersData = vtkSmartPointer::New(); FibersData->SetPoints(vtkSmartPointer::New()); FibersData->SetLines(vtkSmartPointer::New()); FibersData->SetPoints(vNewPoints); FibersData->SetLines(vNewLines); Prediction = mitk::FiberBundle::New(vNewPolyData); // Bundle->SetFiberColors(255, 255, 255); MITK_INFO << "Cells Prediciton"; MITK_INFO << Prediction->GetFiberPolyData()->GetNumberOfCells(); MITK_INFO << "Cells Tractorgram"; MITK_INFO << m_TractogramTest->GetFiberPolyData()->GetNumberOfCells(); return Prediction; } void StreamlineFeatureExtractor::GenerateData() { MITK_INFO << "Update"; // mitk::FiberBundle::Pointer inputPrototypes = mitk::IOUtil::Load("/home/r948e/E132-Projekte/Projects/2022_Peretzke_Interactive_Fiber_Dissection/data/Synt_tract_40_prototypes.trk"); std::vector > T_Prototypes; std::vector > T_TractogramPlus; std::vector > T_TractogramMinus; std::vector > T_TractogramTest; std::vector > T_mergedPrototypes; T_Prototypes = ResampleFibers(m_inputPrototypes); T_TractogramMinus= ResampleFibers(m_TractogramMinus); T_TractogramPlus= ResampleFibers(m_TractogramPlus); /* Merge T_Prototypes, T_TractogramMinus and T_TractogramPlus for extra Features*/ T_mergedPrototypes = MergeTractogram(T_Prototypes, T_TractogramPlus, T_TractogramMinus); MITK_INFO << "Calculate Features"; m_DistancesMinus = CalculateDmdf(T_TractogramMinus, T_mergedPrototypes); m_DistancesPlus = CalculateDmdf(T_TractogramPlus, T_mergedPrototypes); std::ifstream f(m_DistancesTestName); - MITK_INFO << m_DistancesTestName; - if (f.good()) + + if (f.good() && m_activeCycle!=0) { - MITK_INFO << "File exists"; + MITK_INFO << "Loading Features of Tractogram"; m_DistancesTest.clear(); std::ifstream myFile(m_DistancesTestName); if(!myFile.is_open()) throw std::runtime_error("Could not open file"); std::string line; vnl_matrix curline; curline.set_size(1, m_DistancesPlus.at(0).cols()); curline.fill(0.0); float val; while(std::getline(myFile, line)) { // Create a stringstream of the current line std::stringstream ss(line); // MITK_INFO << ss; // Keep track of the current column index int colIdx = 0; // Extract each integer while(ss >> val){ // // Add the current integer to the 'colIdx' column's values vector curline.put(0,colIdx, val); // // If the next token is a comma, ignore it and move on // if(ss.peek() == ',') ss.ignore(); // // Increment the column index colIdx++; } m_DistancesTest.push_back(curline); } // Close file myFile.close(); } else { MITK_INFO << m_DistancesTestName; MITK_INFO << "Resample Test Data"; T_TractogramTest= ResampleFibers(m_TractogramTest); MITK_INFO << "Calculate Features of Test Data"; m_DistancesTest= CalculateDmdf(T_TractogramTest, T_mergedPrototypes); + std::ofstream myFile(m_DistancesTestName); + // myFile << colname << "\n"; + for(long unsigned int i = 0; i < m_DistancesTest.size(); ++i) + { + myFile << m_DistancesTest.at(i); + } + myFile.close(); + } MITK_INFO << m_DistancesTest.size(); MITK_INFO << "Sizes of Plus and Minus"; MITK_INFO << m_DistancesPlus.size() + m_DistancesMinus.size(); + MITK_INFO << "Sizes of Prototypes"; + MITK_INFO << T_mergedPrototypes.size(); + MITK_INFO << T_mergedPrototypes.at(0).rows(); MITK_INFO << "Size of Test Data"; MITK_INFO << m_DistancesTest.size(); MITK_INFO << "Done with Datacreation"; m_index =GetData(); } } - diff --git a/Modules/FiberDissection/MachineLearning/mitkStreamlineFeatureExtractor.h b/Modules/FiberDissection/MachineLearning/mitkStreamlineFeatureExtractor.h index b9dab44..58e23a5 100644 --- a/Modules/FiberDissection/MachineLearning/mitkStreamlineFeatureExtractor.h +++ b/Modules/FiberDissection/MachineLearning/mitkStreamlineFeatureExtractor.h @@ -1,108 +1,110 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #ifndef StreamlineFeatureExtractor_h #define StreamlineFeatureExtractor_h #include "MitkFiberDissectionExports.h" // MITK #include #include // ITK #include // VTK #include #include #include #include #include // OpenCV #include #include #include namespace mitk{ /** * \brief */ class MITKFIBERDISSECTION_EXPORT StreamlineFeatureExtractor { public: StreamlineFeatureExtractor(); ~StreamlineFeatureExtractor(); typedef itk::Image< float, 3 > FloatImageType; typedef itk::Image< unsigned char, 3 > UcharImageType; void Update(){ this->GenerateData(); } void SetTractogramPlus(const mitk::FiberBundle::Pointer &Tractogram); void SetTractogramMinus(const mitk::FiberBundle::Pointer &Tractogram); void SetTractogramTest(const mitk::FiberBundle::Pointer &Tractogram, std::string TractogramTestName); void SetTractogramPrototypes(const mitk::FiberBundle::Pointer &TractogramPrototypes, bool standard); void SetActiveCycle(int &activeCycle); + void SetInitRandom(int &initRandom); // void SetRandomThreshold(int &threshold); void CreateClassifier(); std::vector> GetData(); // void CreatePrediction(std::vector &index); mitk::FiberBundle::Pointer CreatePrediction(std::vector &index); mitk::FiberBundle::Pointer m_Prediction; mitk::FiberBundle::Pointer m_ToLabel; std::vector> m_index; protected: void GenerateData(); std::vector< vnl_matrix > ResampleFibers(FiberBundle::Pointer tractogram); std::vector > CalculateDmdf(std::vector > tractogram, std::vector > prototypes); std::vector< vnl_matrix > MergeTractogram(std::vector > prototypes, std::vector > positive_local_prototypes, std::vector > negative_local_prototypes); std::vector Sort(std::vector sortingVector, int lengths); unsigned int m_NumPoints; int m_activeCycle; + int m_initRandom; int m_thresh; mitk::FiberBundle::Pointer m_TractogramPlus; mitk::FiberBundle::Pointer m_TractogramMinus; mitk::FiberBundle::Pointer m_TractogramTest; mitk::FiberBundle::Pointer m_inputPrototypes; std::string m_DistancesTestName; std::vector > m_DistancesPlus; std::vector > m_DistancesMinus; std::vector > m_DistancesTest; cv::Ptr m_traindata; }; } #endif diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionView.cpp b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionView.cpp index 0bca8ee..aad5e77 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionView.cpp +++ b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionView.cpp @@ -1,1196 +1,1233 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ // Blueberry #include #include #include // Qmitk #include "QmitkInteractiveFiberDissectionView.h" #include //Pointset #include //Pointset #include #include #include #include #include #include #include #include "mitkNodePredicateDataType.h" #include #include #include #include //#include #include #include "usModuleRegistry.h" //#include #include #include #include #include #include #include #include #include #include #include #include #include #include const std::string QmitkInteractiveFiberDissectionView::VIEW_ID = "org.mitk.views.interactivefiberdissection"; const std::string id_DataManager = "org.mitk.views.datamanager"; using namespace mitk; QmitkInteractiveFiberDissectionView::QmitkInteractiveFiberDissectionView() : QmitkAbstractView() , m_Controls( 0 ) , m_IterationCounter(0) , m_RandomExtractionCounter(0) , m_activeCycleCounter(0) , m_StreamlineInteractor(nullptr) { } // Destructor QmitkInteractiveFiberDissectionView::~QmitkInteractiveFiberDissectionView() { //disable interactor if (m_StreamlineInteractor != nullptr) { // m_StreamlineInteractor->SetStreamlineNode(nullptr); m_StreamlineInteractor->EnableInteraction(false); } } void QmitkInteractiveFiberDissectionView::CreateQtPartControl( QWidget *parent ) { // build up qt view, unless already done if ( !m_Controls ) { // create GUI widgets from the Qt Designer's .ui file m_Controls = new Ui::QmitkInteractiveFiberDissectionViewControls; m_Controls->setupUi( parent ); m_Controls->m_selectedPointSetWidget->SetDataStorage(GetDataStorage());//pointset m_Controls->m_selectedPointSetWidget->SetNodePredicate(mitk::NodePredicateAnd::New(//pointset mitk::TNodePredicateDataType::New(),//pointset mitk::NodePredicateNot::New(mitk::NodePredicateOr::New(//pointset mitk::NodePredicateProperty::New("helper object"),//pointset mitk::NodePredicateProperty::New("hidden object")))));//pointset m_Controls->m_selectedPointSetWidget->SetSelectionIsOptional(true);//pointset m_Controls->m_selectedPointSetWidget->SetAutoSelectNewNodes(true);//pointset m_Controls->m_selectedPointSetWidget->SetEmptyInfo(QString("Please select a point set"));//pointset m_Controls->m_selectedPointSetWidget->SetPopUpTitel(QString("Select point set"));//pointsett m_Controls->m_BundleBox->SetDataStorage(this->GetDataStorage()); mitk::TNodePredicateDataType::Pointer isBundle= mitk::TNodePredicateDataType::New(); m_Controls->m_BundleBox->SetPredicate( isBundle ); m_Controls->m_PrototypeBox->SetDataStorage(this->GetDataStorage()); mitk::TNodePredicateDataType::Pointer isPrototype = mitk::TNodePredicateDataType::New(); m_Controls->m_PrototypeBox->SetPredicate( isPrototype ); connect(m_Controls->m_ErazorButton, SIGNAL(toggled(bool)), this, SLOT( RemovefromBundle(bool) ) ); //need connect(m_Controls->m_StreamlineCreation, SIGNAL( clicked() ), this, SLOT( CreateStreamline())); connect(m_Controls->m_AddRandomFibers, SIGNAL( clicked() ), this, SLOT( ExtractRandomFibersFromTractogram() ) ); //need connect(m_Controls->m_TrainClassifier, SIGNAL( clicked() ), this, SLOT( StartAlgorithm( ))); connect(m_Controls->m_CreatePrediction, SIGNAL( clicked() ), this, SLOT( CreatePredictionNode( ))); connect(m_Controls->m_AddUncertainFibers, SIGNAL( clicked() ), this, SLOT( CreateUncertaintySampleNode( ))); connect(m_Controls->m_AddDistanceFibers, SIGNAL( clicked() ), this, SLOT( CreateDistanceSampleNode( ))); connect(m_Controls->m_unclabeling, SIGNAL(toggled(bool)), this, SLOT( RemovefromUncertainty(bool) ) ); //need connect(m_Controls->m_distlabeling, SIGNAL(toggled(bool)), this, SLOT( RemovefromDistance(bool) ) ); //need connect(m_Controls->m_predlabeling, SIGNAL(toggled(bool)), this, SLOT( RemovefromPrediction(bool) ) ); //need + + connect(m_Controls->m_sellabeling, SIGNAL(toggled(bool)), this, SLOT( RemovefromSelection(bool) ) ); //need connect(m_Controls->m_ResampleButton, SIGNAL( clicked() ), this, SLOT( ResampleTractogram( ) ) ); connect(m_Controls->m_RandomPrototypesButton, SIGNAL( clicked() ), this, SLOT( RandomPrototypes( ) ) ); connect(m_Controls->m_SFFPrototypesButton, SIGNAL( clicked() ), this, SLOT( SFFPrototypes( ) ) ); connect(m_Controls->m_addPointSetPushButton, &QPushButton::clicked,//pointset this, &QmitkInteractiveFiberDissectionView::OnAddPointSetClicked);//pointset connect(m_Controls->m_selectedPointSetWidget, &QmitkSingleNodeSelectionWidget::CurrentSelectionChanged,//pointset this, &QmitkInteractiveFiberDissectionView::OnCurrentSelectionChanged);//pointset auto renderWindowPart = this->GetRenderWindowPart();//pointset if (nullptr != renderWindowPart)//pointset this->RenderWindowPartActivated(renderWindowPart);//pointset this->OnCurrentSelectionChanged(m_Controls->m_selectedPointSetWidget->GetSelectedNodes());//pointset } UpdateGui(); } void QmitkInteractiveFiberDissectionView::SetFocus() { m_Controls->toolBoxx->setFocus(); //m_Controls->m_addPointSetPushButton->setFocus();//pointset } void QmitkInteractiveFiberDissectionView::UpdateGui() { m_Controls->m_FibLabel->setText("mandatory"); m_Controls->m_InputData->setTitle("Please Select Input Data"); // disable alle frames m_Controls->m_ErazorButton->setCheckable(true); m_Controls->m_ErazorButton->setEnabled(false); m_Controls->m_unclabeling->setCheckable(true); m_Controls->m_unclabeling->setEnabled(false); m_Controls->m_predlabeling->setCheckable(true); m_Controls->m_predlabeling->setEnabled(false); m_Controls->m_distlabeling->setCheckable(true); m_Controls->m_distlabeling->setEnabled(false); + m_Controls->m_sellabeling->setCheckable(true); + m_Controls->m_sellabeling->setEnabled(false); m_Controls->m_addPointSetPushButton->setEnabled(false); m_Controls->m_StreamlineCreation->setEnabled(false); m_Controls->m_TrainClassifier->setEnabled(false); m_Controls->m_CreatePrediction->setEnabled(false); m_Controls->m_CreateUncertantyMap->setEnabled(false); m_Controls->m_Numtolabel->setEnabled(false); m_Controls->m_Numtolabel2->setEnabled(false); m_Controls->m_addPointSetPushButton->setEnabled(false); m_Controls->m_AddRandomFibers->setEnabled(false); m_Controls->m_AddDistanceFibers->setEnabled(false); m_Controls->m_AddUncertainFibers->setEnabled(false); - m_Controls->m_unclabeling->setEnabled(false); - m_Controls->m_predlabeling->setEnabled(false); - m_Controls->m_distlabeling->setEnabled(false); // m_Controls->m_PrototypeBox->setEditable(false); // m_Controls->m_useStandardP-> bool fibSelected = !m_SelectedFB.empty(); bool multipleFibsSelected = (m_SelectedFB.size()>1); bool sthSelected = m_SelectedImageNode.IsNotNull(); bool psSelected = m_SelectedPS.IsNotNull(); // bool nfibSelected = !m_newfibersBundleNode.empty(); // bool posSelected = !m_positivBundlesNode.empty(); bool nfibSelected = m_newfibersBundleNode.IsNotNull(); // bool posSelected = !m_positivBundlesNode.IsNotNull(); // bool negSelected = !m_negativeBundleNode.IsNotNull(); bool posSelected = this->GetDataStorage()->Exists(m_positivBundlesNode); bool negSelected = this->GetDataStorage()->Exists(m_negativeBundleNode); bool indexSelected = !m_index.empty(); bool uncertaintySelected = this->GetDataStorage()->Exists(m_UncertaintyLabelNode); bool distanceSelected = this->GetDataStorage()->Exists(m_DistanceLabelNode); bool predictionSelected = this->GetDataStorage()->Exists(m_PredictionNode); // toggle visibility of elements according to selected method // are fiber bundles selected? if ( fibSelected ) { m_Controls->m_FibLabel->setText(QString(m_SelectedFB.at(0)->GetName().c_str())); m_Controls->m_addPointSetPushButton->setEnabled(true); m_Controls->m_AddRandomFibers->setEnabled(true); + m_Controls->m_sellabeling->setEnabled(true); // more than two bundles needed to join/subtract if (multipleFibsSelected) { m_Controls->m_FibLabel->setText("multiple bundles selected"); } } // is image selected if (sthSelected) { m_Controls->m_addPointSetPushButton->setEnabled(true); } if (psSelected) { m_Controls->m_StreamlineCreation->setEnabled(true); } if (nfibSelected && posSelected) { m_Controls->m_ErazorButton->setEnabled(true); } if (posSelected && negSelected) { m_Controls->m_TrainClassifier->setEnabled(true); } if (indexSelected) { m_Controls->m_CreatePrediction->setEnabled(true); m_Controls->m_AddUncertainFibers->setEnabled(true); m_Controls->m_Numtolabel->setEnabled(true); m_Controls->m_AddDistanceFibers->setEnabled(true); m_Controls->m_Numtolabel2->setEnabled(true); } if (uncertaintySelected) { m_Controls->m_unclabeling->setEnabled(true); } if (predictionSelected) { m_Controls->m_predlabeling->setEnabled(true); } if (distanceSelected) { m_Controls->m_distlabeling->setEnabled(true); } - if (m_Controls->m_useStandardP->isChecked()) - { - m_Controls->m_PrototypeBox->setEditable(true); - } +// if (m_Controls->m_useStandardP->isChecked()) +// { +// m_Controls->m_PrototypeBox->setEditable(true); +// } } void QmitkInteractiveFiberDissectionView::OnEndInteraction() { } void QmitkInteractiveFiberDissectionView::ResampleTractogram() { mitk::DataNode::Pointer node = m_Controls->m_BundleBox->GetSelectedNode(); auto tractogram = dynamic_cast(node->GetData()); mitk::FiberBundle::Pointer temp_fib = tractogram->GetDeepCopy(); temp_fib->ResampleToNumPoints(40); MITK_INFO << "Resampling Done"; mitk::DataNode::Pointer newnode = mitk::DataNode::New(); newnode->SetData( temp_fib ); newnode->SetName( node->GetName() + "_" + std::to_string(40) ); this->GetDataStorage()->Add(newnode); UpdateGui(); } void QmitkInteractiveFiberDissectionView::RandomPrototypes() { MITK_INFO << "Number of Fibers to use as Prototypes: "; MITK_INFO << m_Controls->m_NumPrototypes->value(); mitk::FiberBundle::Pointer fib = dynamic_cast(m_Controls->m_BundleBox->GetSelectedNode()->GetData()); MITK_INFO << fib->GetNumFibers(); std::vector myvec; for (unsigned int k=0; kGetNumFibers(); k++) { myvec.push_back(k); } // auto rng = std::default_random_engine {}; -// std::random_shuffle(std::begin(myvec), std::end(myvec)); + std::random_shuffle(std::begin(myvec), std::end(myvec)); vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); vtkSmartPointer weights = vtkSmartPointer::New(); /* Check wether all Streamlines of the bundles are labeled... If all are labeled Skip for Loop*/ unsigned int counter = 0; for (int i=0; im_NumPrototypes->value(); i++) { vtkCell* cell = fib->GetFiberPolyData()->GetCell(myvec.at(i)); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights->InsertValue(counter, fib->GetFiberWeight(myvec.at(i))); vNewLines->InsertNextCell(container); counter++; } vNewPolyData->SetLines(vNewLines); vNewPolyData->SetPoints(vNewPoints); mitk::FiberBundle::Pointer PrototypesBundle = mitk::FiberBundle::New(vNewPolyData); PrototypesBundle->SetFiberWeights(weights); mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetData(PrototypesBundle); node->SetName("Random_Prototypes"); // MITK_INFO << "Number of Streamlines in first function"; // MITK_INFO << m_newfibersBundleNode->GetData()->GetFiberPolyData()->GetNumberOfCells(); this->GetDataStorage()->Add(node); } void QmitkInteractiveFiberDissectionView::SFFPrototypes() { MITK_INFO << "Number of Fibers to use as Prototypes: "; MITK_INFO << m_Controls->m_NumPrototypes->value(); MITK_INFO << "Start Creating Prototypes based on SFF"; mitk::FiberBundle::Pointer fib = dynamic_cast(m_Controls->m_BundleBox->GetSelectedNode()->GetData()); /* Get Subset of Tractogram*/ int size_subset = std::max(1.0, ceil(3.0 * m_Controls->m_NumPrototypes->value() * std::log(m_Controls->m_NumPrototypes->value()))); MITK_INFO << fib->GetNumFibers(); std::vector myvec; for (unsigned int k=0; kGetNumFibers(); k++) { myvec.push_back(k); } // std::random_shuffle(std::begin(myvec), std::end(myvec)); vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); vtkSmartPointer weights = vtkSmartPointer::New(); unsigned int counter = 0; for (int i=0; iGetFiberPolyData()->GetCell(myvec.at(i)); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights->InsertValue(counter, fib->GetFiberWeight(myvec.at(i))); vNewLines->InsertNextCell(container); counter++; } vNewPolyData->SetLines(vNewLines); vNewPolyData->SetPoints(vNewPoints); mitk::FiberBundle::Pointer temp_fib = mitk::FiberBundle::New(vNewPolyData); temp_fib->SetFiberWeights(weights); MITK_INFO << temp_fib->GetFiberPolyData()->GetNumberOfCells(); /* Create std::vector of the SubsetBundle*/ std::vector< vnl_matrix > out_fib(temp_fib->GetFiberPolyData()->GetNumberOfCells()); for (int i=0; iGetFiberPolyData()->GetNumberOfCells(); i++) { vtkCell* cell = temp_fib->GetFiberPolyData()->GetCell(i); int numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vnl_matrix streamline; streamline.set_size(3, cell->GetNumberOfPoints()); streamline.fill(0.0); for (int j=0; jGetPoint(j, cand); vnl_vector_fixed< float, 3 > candV; candV[0]=cand[0]; candV[1]=cand[1]; candV[2]=cand[2]; streamline.set_column(j, candV); } // out_fib.push_back(streamline); out_fib.at(i)=streamline; } /* Calculate the distancematrix of Subset*/ std::vector< vnl_matrix > dist_vec(out_fib.size());// cv::parallel_for_(cv::Range(0, out_fib.size()), [&](const cv::Range &range) { for (int i = range.start; i < range.end; i++) // for (unsigned int i=0; i distances; distances.set_size(1, out_fib.size()); distances.fill(0.0); for (unsigned int j=0; j single_distances; single_distances.set_size(1, out_fib.at(0).cols()); single_distances.fill(0.0); vnl_matrix single_distances_flip; single_distances_flip.set_size(1, out_fib.at(0).cols()); single_distances_flip.fill(0.0); for (unsigned int ik=0; ik single_distances.mean()) { distances.put(0,j, single_distances.mean()); } else { distances.put(0,j, single_distances_flip.mean()); } } // dist_vec.push_back(distances); dist_vec.at(i) = distances; } }); /*Index to find values in distancematrix*/ std::vector myidx; /*Index to find actual streamlines using indexUnc*/ std::vector indexUncDist; /*Start with the Streamline of the highest entropy, which is in distance_matrix at idx 0*/ myidx.push_back(0); /*Vecotr that stores minvalues of current iteration*/ vnl_matrix cur_vec; cur_vec.set_size(1, size_subset); cur_vec.fill(0.0); for (int i=0; im_NumPrototypes->value(); i++) { // unsigned int cur_i = indexUnc.at(myidx.at(i)); /*Save mean distance of all used Samples*/ vnl_matrix sum_matrix; sum_matrix.set_size(myidx.size(), size_subset); sum_matrix.fill(0); for (unsigned int ii=0; ii vNewPolyData2 = vtkSmartPointer::New(); vtkSmartPointer vNewLines2 = vtkSmartPointer::New(); vtkSmartPointer vNewPoints2 = vtkSmartPointer::New(); vtkSmartPointer weights2 = vtkSmartPointer::New(); /* Check wether all Streamlines of the bundles are labeled... If all are labeled Skip for Loop*/ counter = 0; for (int i=0; im_NumPrototypes->value(); i++) { vtkCell* cell = fib->GetFiberPolyData()->GetCell(myidx.at(i)); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints2->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights2->InsertValue(counter, fib->GetFiberWeight(myvec.at(i))); vNewLines2->InsertNextCell(container); counter++; } vNewPolyData2->SetLines(vNewLines2); vNewPolyData2->SetPoints(vNewPoints2); mitk::FiberBundle::Pointer PrototypesBundle = mitk::FiberBundle::New(vNewPolyData2); PrototypesBundle->SetFiberWeights(weights2); MITK_INFO << PrototypesBundle->GetFiberPolyData()->GetNumberOfCells(); mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetData(PrototypesBundle); node->SetName("SFF_Prototypes"); //// MITK_INFO << "Number of Streamlines in first function"; //// MITK_INFO << m_newfibersBundleNode->GetData()->GetFiberPolyData()->GetNumberOfCells(); this->GetDataStorage()->Add(node); } void QmitkInteractiveFiberDissectionView::OnAddPointSetClicked()//pointset { // ask for the name of the point set bool ok = false; QString name = QInputDialog::getText(QApplication::activeWindow(), tr("Add point set..."), tr("Enter name for the new point set"), QLineEdit::Normal, tr("PointSet").arg(++m_IterationCounter), &ok); // QString name = "PointSet"; if (!ok || name.isEmpty()) { return; } mitk::PointSet::Pointer pointSet = mitk::PointSet::New(); mitk::DataNode::Pointer pointSetNode = mitk::DataNode::New(); pointSetNode->SetData(pointSet); pointSetNode->SetProperty("name", mitk::StringProperty::New(name.toStdString())); pointSetNode->SetProperty("opacity", mitk::FloatProperty::New(1)); pointSetNode->SetColor(1.0, 1.0, 0.0); this->GetDataStorage()->Add(pointSetNode, m_SelectedImageNode); m_Controls->m_selectedPointSetWidget->SetCurrentSelectedNode(pointSetNode); } void QmitkInteractiveFiberDissectionView::OnCurrentSelectionChanged(QmitkSingleNodeSelectionWidget::NodeList /*nodes*/)//pointset { m_Controls->m_poinSetListWidget->SetPointSetNode(m_Controls->m_selectedPointSetWidget->GetSelectedNode()); m_SelectedPS = m_Controls->m_selectedPointSetWidget->GetSelectedNode(); // m_Controls->m_trainbundleWidget->SetPointSetNode(m_Controls->m_trainbundleWidget->GetSelectedNode()); // m_trainbundle = m_Controls->m_trainbundleWidget->GetSelectedNode(); UpdateGui(); } void QmitkInteractiveFiberDissectionView::OnSelectionChanged(berry::IWorkbenchPart::Pointer /*part*/, const QList& nodes) { m_SelectedFB.clear(); if (nodes.empty() || nodes.front().IsNull()) { m_SelectedImageNode = nullptr; } else { m_SelectedImageNode = nodes.front(); } for (auto node: nodes) { if (dynamic_cast(node->GetData())) m_SelectedImage = dynamic_cast(node->GetData()); else if ( dynamic_cast(node->GetData()) ) m_SelectedFB.push_back(node); } UpdateGui(); } void QmitkInteractiveFiberDissectionView::RenderWindowPartActivated(mitk::IRenderWindowPart* renderWindowPart)//pointset { if (nullptr != m_Controls) { m_Controls->m_poinSetListWidget->AddSliceNavigationController(renderWindowPart->GetQmitkRenderWindow("axial")->GetSliceNavigationController()); m_Controls->m_poinSetListWidget->AddSliceNavigationController(renderWindowPart->GetQmitkRenderWindow("sagittal")->GetSliceNavigationController()); m_Controls->m_poinSetListWidget->AddSliceNavigationController(renderWindowPart->GetQmitkRenderWindow("coronal")->GetSliceNavigationController()); } } void QmitkInteractiveFiberDissectionView::RenderWindowPartDeactivated(mitk::IRenderWindowPart* renderWindowPart)//pointset { if (nullptr != m_Controls) { m_Controls->m_poinSetListWidget->RemoveSliceNavigationController(renderWindowPart->GetQmitkRenderWindow("axial")->GetSliceNavigationController()); m_Controls->m_poinSetListWidget->RemoveSliceNavigationController(renderWindowPart->GetQmitkRenderWindow("sagittal")->GetSliceNavigationController()); m_Controls->m_poinSetListWidget->RemoveSliceNavigationController(renderWindowPart->GetQmitkRenderWindow("coronal")->GetSliceNavigationController()); } } void QmitkInteractiveFiberDissectionView::CreateStreamline() { if (m_positivBundlesNode.IsNull()) { mitk::DataNode::Pointer node = mitk::DataNode::New(); m_positiveFibersData = vtkSmartPointer::New(); m_positiveFibersData->SetPoints(vtkSmartPointer::New()); m_positiveFibersData->SetLines(vtkSmartPointer::New()); m_positiveBundle = mitk::FiberBundle:: New(m_positiveFibersData); node->SetData( m_positiveBundle ); m_positivBundlesNode = node; this->GetDataStorage()->Add(m_positivBundlesNode); MITK_INFO << "Create Bundle"; } if (!m_positivBundlesNode.IsNull()) { this->GetDataStorage()->Remove(m_positivBundlesNode); MITK_INFO << "Adding fibers"; MITK_INFO << m_positiveBundle->GetFiberPolyData()->GetNumberOfCells(); m_positiveFibersData = m_positiveBundle->GetFiberPolyData(); } vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); unsigned int counter = 0; for (unsigned int i=0; iGetNumberOfCells(); ++i) { MITK_INFO<< "New Line"; vtkCell* cell = m_positiveFibersData->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } vNewLines->InsertNextCell(container); counter++; } mitk::PointSet::Pointer pointSet = dynamic_cast(m_SelectedPS->GetData()); vnl_matrix streamline; streamline.set_size(3, pointSet->GetSize()); streamline.fill(0.0); mitk::PointSet::PointsIterator begin = pointSet->Begin(); mitk::PointSet::PointsIterator end = pointSet->End(); unsigned int i; mitk::PointSet::PointsContainer::Iterator it; for (it = begin, i = 0; it != end; ++it, ++i) { PointSet::PointType pt = pointSet->GetPoint(it->Index()); vnl_vector_fixed< float, 3 > candV; candV[0]=pt[0]; candV[1]=pt[1]; candV[2]=pt[2]; streamline.set_column(i, candV); } // build Fiber vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jInsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } MITK_INFO<< "Last Line from current pointset"; vNewLines->InsertNextCell(container); vNewPolyData->SetPoints(vNewPoints); vNewPolyData->SetLines(vNewLines); m_positiveFibersData = vtkSmartPointer::New(); m_positiveFibersData->SetPoints(vtkSmartPointer::New()); m_positiveFibersData->SetLines(vtkSmartPointer::New()); m_positiveFibersData->SetPoints(vNewPoints); m_positiveFibersData->SetLines(vNewLines); m_positiveBundle = mitk::FiberBundle::New(vNewPolyData); // m_positiveBundle->SetTrackVisHeader(dynamic_cast(m_SelectedImageNode->GetData())->GetGeometry()); m_positiveBundle->SetFiberColors(0, 255, 0); mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetData(m_positiveBundle); node->SetName("+Bundle"); m_positivBundlesNode= node; MITK_INFO << "The + Bundle has Streamlines:"; auto m_PosStreamline= dynamic_cast(m_positivBundlesNode->GetData()); MITK_INFO << m_PosStreamline->GetFiberPolyData()->GetNumberOfCells(); this->GetDataStorage()->Add(m_positivBundlesNode); // m_Controls->m_selectedPointSetWidget->m_ToggleAddPoint->setEnabled(false); UpdateGui(); } void QmitkInteractiveFiberDissectionView::ExtractRandomFibersFromTractogram() { m_SelectedFB.at(0)->SetVisibility(false); m_Controls->m_ErazorButton->setChecked(false); MITK_INFO << "Number of Fibers to extract from Tractogram: "; MITK_INFO << m_Controls->m_NumRandomFibers->value(); if (this->GetDataStorage()->Exists(m_newfibersBundleNode)) { MITK_INFO << "To Label Bundle Exists"; mitk::FiberBundle::Pointer Stack = dynamic_cast(m_newfibersBundleNode->GetData()); this->GetDataStorage()->Remove(m_newfibersBundleNode); mitk::DataNode::Pointer node = mitk::DataNode::New(); m_newfibersFibersData = vtkSmartPointer::New(); m_newfibersFibersData->SetPoints(vtkSmartPointer::New()); m_newfibersBundle = mitk::FiberBundle:: New(m_newfibersFibersData); m_newfibersFibersData->SetLines(vtkSmartPointer::New()); // node->SetData( m_newfibersBundle ); // m_newfibersBundleNode = node ; MITK_INFO << "Create Bundle"; } mitk::FiberBundle::Pointer fib = dynamic_cast(m_SelectedFB.at(0)->GetData()); // mitk::FiberBundle::Pointer fib = dynamic_cast(m_trainbundle->GetData()); vtkSmartPointer vNewPolyData = vtkSmartPointer::New(); vtkSmartPointer vNewLines = vtkSmartPointer::New(); vtkSmartPointer vNewPoints = vtkSmartPointer::New(); vtkSmartPointer weights = vtkSmartPointer::New(); // weights->SetNumberOfValues(this->GetNumFibers()+fib->GetNumFibers()); /* Check weather all Streamlines of the bundles are labeled... If all are labeled Skip for Loop*/ unsigned int counter = 0; int thresh1; int thresh2; thresh2 = m_Controls->m_NumRandomFibers->value()*(m_RandomExtractionCounter+1); thresh1 = m_Controls->m_NumRandomFibers->value()*(m_RandomExtractionCounter); if (thresh1>fib->GetFiberPolyData()->GetNumberOfCells()) { thresh1=fib->GetFiberPolyData()->GetNumberOfCells(); } if (thresh2>fib->GetFiberPolyData()->GetNumberOfCells()) { thresh2=fib->GetFiberPolyData()->GetNumberOfCells(); } if (thresh1!=fib->GetFiberPolyData()->GetNumberOfCells()) { for ( int i=thresh1; iGetFiberPolyData()->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vtkSmartPointer container = vtkSmartPointer::New(); for (unsigned int j=0; jGetPoint(j, p); vtkIdType id = vNewPoints->InsertNextPoint(p); container->GetPointIds()->InsertNextId(id); } weights->InsertValue(counter, fib->GetFiberWeight(i)); vNewLines->InsertNextCell(container); counter++; } vNewPolyData->SetLines(vNewLines); vNewPolyData->SetPoints(vNewPoints); m_newfibersFibersData = vtkSmartPointer::New(); m_newfibersFibersData->SetPoints(vtkSmartPointer::New()); m_newfibersFibersData->SetLines(vtkSmartPointer::New()); m_newfibersFibersData->SetPoints(vNewPoints); m_newfibersFibersData->SetLines(vNewLines); m_newfibersBundle = mitk::FiberBundle::New(vNewPolyData); m_newfibersBundle->SetFiberColors(255, 255, 255); m_newfibersBundle->SetFiberWeights(weights); mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetData(m_newfibersBundle); node->SetName("ToLabel"); m_newfibersBundleNode = node; // MITK_INFO << "Number of Streamlines in first function"; // MITK_INFO << m_newfibersBundleNode->GetData()->GetFiberPolyData()->GetNumberOfCells(); this->GetDataStorage()->Add(m_newfibersBundleNode); m_RandomExtractionCounter++; } m_Controls->m_ErazorButton->setChecked(true); UpdateGui(); } void QmitkInteractiveFiberDissectionView::RemovefromBundle( bool checked ) { if (checked) { if (m_StreamlineInteractor.IsNull()) { this->CreateStreamlineInteractor(); mitk::FiberBundle::Pointer m_negativeBundle = mitk::FiberBundle::New(); mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetName("-Bundle"); node->SetData(m_negativeBundle); m_negativeBundleNode = node; this->GetDataStorage()->Add(m_negativeBundleNode); m_StreamlineInteractor->EnableInteraction(true); + m_StreamlineInteractor->LabelfromPrediction(false); m_StreamlineInteractor->SetNegativeNode(m_negativeBundleNode); m_StreamlineInteractor->SetPositiveNode(m_positivBundlesNode); m_StreamlineInteractor->SetToLabelNode(m_newfibersBundleNode); } else { m_StreamlineInteractor->EnableInteraction(true); + m_StreamlineInteractor->LabelfromPrediction(false); m_StreamlineInteractor->SetPositiveNode(m_positivBundlesNode); m_StreamlineInteractor->SetToLabelNode(m_newfibersBundleNode); } } else { m_StreamlineInteractor->EnableInteraction(false); // m_StreamlineInteractor = nullptr; } UpdateGui(); } void QmitkInteractiveFiberDissectionView::CreateStreamlineInteractor() { m_StreamlineInteractor = mitk::StreamlineInteractor::New(); m_StreamlineInteractor->LoadStateMachine("Streamline3DStates.xml", us::ModuleRegistry::GetModule("MitkFiberDissection")); m_StreamlineInteractor->SetEventConfig("Streamline3DConfig.xml", us::ModuleRegistry::GetModule("MitkFiberDissection")); // m_StreamlineInteractor->SetRotationEnabled(rotationEnabled); } void QmitkInteractiveFiberDissectionView::StartAlgorithm() { + m_negativeBundle = dynamic_cast(m_negativeBundleNode->GetData()); + if (m_activeCycleCounter==0) + { + + int neg_cells; + int pos_cells; + neg_cells = m_negativeBundle->GetFiberPolyData()->GetNumberOfCells(); + pos_cells = m_positiveBundle->GetFiberPolyData()->GetNumberOfCells(); + m_initRandom = neg_cells - pos_cells; + MITK_INFO << m_initRandom; + } + this->GetDataStorage()->Remove(m_UncertaintyLabelNode); this->GetDataStorage()->Remove(m_DistanceLabelNode); // vtkSmartPointer< vtkFloatArray > weights = m_positiveBundle->GetFiberWeights(); // MITK_INFO << weights->GetNumberOfValues(); // MITK_INFO << m_positiveBundle->GetNumFibers(); // for (unsigned int k=0; kGetNumFibers(); k++) // { // m_positiveBundle->SetFiberWeight(k, 1.0); // } // vtkSmartPointer< vtkFloatArray > newweights = m_positiveBundle->GetFiberWeights(); // MITK_INFO << newweights->GetNumberOfValues(); m_Controls->m_unclabeling->setChecked(false); m_Controls->m_distlabeling->setChecked(false); m_Controls->m_predlabeling->setChecked(false); clusterer.reset(); MITK_INFO << "Extract Features"; - m_negativeBundle = dynamic_cast(m_negativeBundleNode->GetData()); + clusterer = std::make_shared(); + clusterer->SetInitRandom(m_initRandom); clusterer->SetActiveCycle(m_activeCycleCounter); clusterer->SetTractogramPlus(m_positiveBundle); clusterer->SetTractogramMinus(m_negativeBundle); clusterer->SetTractogramPrototypes(dynamic_cast(m_Controls->m_PrototypeBox->GetSelectedNode()->GetData()), m_Controls->m_useStandardP->isChecked()); clusterer->SetTractogramTest(dynamic_cast(m_SelectedFB.at(0)->GetData()), m_SelectedFB.at(0)->GetName()); // clusterer->SetTractogramTest(dynamic_cast(m_trainbundle->GetData()), m_trainbundle->GetName()); clusterer->Update(); m_index = clusterer->m_index; MITK_INFO << "Number of Cycles"; MITK_INFO << m_activeCycleCounter; m_activeCycleCounter += 1; // m_Prediction = clusterer->CreatePrediction(m_index.at(0)); // mitk::DataNode::Pointer node = mitk::DataNode::New(); // node->SetData(m_Prediction); // node->SetName("Prediction"); // m_PredictionNode = node; // this->GetDataStorage()->Add(m_PredictionNode); // m_UncertaintyLabel = clusterer->m_UncertaintyLabel; // mitk::DataNode::Pointer node2 = mitk::DataNode::New(); // node2->SetData(m_UncertaintyLabel); // node2->SetName("UncertaintyLabels"); // m_UncertaintyLabelNode = node2; // MITK_INFO << "Number of Streamlines in first function"; // MITK_INFO << m_newfibersBundleNode->GetData()->GetFiberPolyData()->GetNumberOfCells(); // this->GetDataStorage()->Add(m_UncertaintyLabelNode); // this->GetDataStorage()->Add(m_PredictionNode); // clusterer->GetData(); // MITK_INFO << data.at(0); // MITK_INFO << data.at(1); // cv::Ptr m_traindata = clusterer->GetData(); // MITK_INFO << clusterer->m_labels; // MITK_INFO << data.at(1); // MITK_INFO << "Start Classification"; // clusterer->CreateClassifier(); // cv::Mat curdata = clusterer->StartAlgorithm(); // MITK_INFO << curdata; MITK_INFO << "Algorithm run succesfully"; m_Controls->m_CreatePrediction->setEnabled(true); UpdateGui(); } void QmitkInteractiveFiberDissectionView::CreatePredictionNode() { MITK_INFO << "Create Prediction"; m_Prediction = clusterer->CreatePrediction(m_index.at(0)); mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetData(m_Prediction); auto s = std::to_string(m_activeCycleCounter); node->SetName("Prediction"+s); m_PredictionNode = node; this->GetDataStorage()->Add(m_PredictionNode); UpdateGui(); } void QmitkInteractiveFiberDissectionView::CreateUncertaintySampleNode() { MITK_INFO << "Create Fibers to label based on Uncertainty"; std::vector myvec = m_index.at(1); myvec.resize(m_Controls->m_Numtolabel->value()); MITK_INFO << m_index.at(1).size(); MITK_INFO << myvec.size(); m_UncertaintyLabel = clusterer->CreatePrediction(myvec); mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetData(m_UncertaintyLabel); auto s = std::to_string(m_activeCycleCounter); node->SetName("UncertaintyLabel"+s); m_UncertaintyLabelNode = node; this->GetDataStorage()->Add(m_UncertaintyLabelNode); UpdateGui(); } void QmitkInteractiveFiberDissectionView::CreateDistanceSampleNode() { MITK_INFO << "Create Fibers to label based on Distance in Features-Space"; std::vector myvec = m_index.at(2); myvec.resize(m_Controls->m_Numtolabel2->value()); MITK_INFO << m_index.at(2).size(); MITK_INFO << myvec.size(); m_DistanceLabel = clusterer->CreatePrediction(myvec); mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetData(m_DistanceLabel); auto s = std::to_string(m_activeCycleCounter); node->SetName("DistanceLabel"+s); m_DistanceLabelNode = node; this->GetDataStorage()->Add(m_DistanceLabelNode); UpdateGui(); } void QmitkInteractiveFiberDissectionView::RemovefromUncertainty( bool checked ) { if (checked) { m_UncertaintyLabel->SetFiberColors(255, 255, 255); m_StreamlineInteractor->EnableInteraction(true); + m_StreamlineInteractor->LabelfromPrediction(false); m_StreamlineInteractor->SetToLabelNode(m_UncertaintyLabelNode); } else { m_StreamlineInteractor->EnableInteraction(false); // m_StreamlineInteractor = nullptr; } RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkInteractiveFiberDissectionView::RemovefromDistance( bool checked ) { if (checked) { m_DistanceLabel->SetFiberColors(255, 255, 255); m_StreamlineInteractor->EnableInteraction(true); + m_StreamlineInteractor->LabelfromPrediction(false); m_StreamlineInteractor->SetToLabelNode(m_DistanceLabelNode); } else { m_StreamlineInteractor->EnableInteraction(false); // m_StreamlineInteractor = nullptr; } RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkInteractiveFiberDissectionView::RemovefromPrediction( bool checked ) { if (checked) { // m_Prediction->SetFiberColors(255, 255, 255); m_StreamlineInteractor->EnableInteraction(true); + m_StreamlineInteractor->LabelfromPrediction(true); m_StreamlineInteractor->SetToLabelNode(m_PredictionNode); } else { m_StreamlineInteractor->EnableInteraction(false); // m_StreamlineInteractor = nullptr; } } + +void QmitkInteractiveFiberDissectionView::RemovefromSelection( bool checked ) +{ + if (checked) + { + + // m_Prediction->SetFiberColors(255, 255, 255); + m_StreamlineInteractor->EnableInteraction(true); + m_StreamlineInteractor->LabelfromPrediction(true); + m_StreamlineInteractor->SetToLabelNode(m_SelectedFB.at(0)); + } + else + { + m_StreamlineInteractor->EnableInteraction(false); + // m_StreamlineInteractor = nullptr; + } +} diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionView.h b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionView.h index 1712ea2..7573d57 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionView.h +++ b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionView.h @@ -1,161 +1,163 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #ifndef QmitkInteractiveFiberDissectionView_h #define QmitkInteractiveFiberDissectionView_h #include "ui_QmitkInteractiveFiberDissectionViewControls.h" #include //Pointset #include //Pointset #include //Pointset #include //Pointset #include #include #include //Pointset #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include /*! \brief View to process fiber bundles. Supplies methods to extract fibers from the bundle, fiber resampling, mirroring, join and subtract bundles and much more. */ class QmitkInteractiveFiberDissectionView : public QmitkAbstractView, public mitk::IRenderWindowPartListener { // this is needed for all Qt objects that should have a Qt meta-object // (everything that derives from QObject and wants to have signal/slots) Q_OBJECT public: typedef itk::Image< unsigned char, 3 > ItkUCharImageType; typedef itk::Image< float, 3 > ItkFloatImageType; static const std::string VIEW_ID; QmitkInteractiveFiberDissectionView(); virtual ~QmitkInteractiveFiberDissectionView(); virtual void CreateQtPartControl(QWidget *parent) override; /// /// Sets the focus to an internal widget. /// virtual void SetFocus() override; protected slots: void RenderWindowPartActivated(mitk::IRenderWindowPart* renderWindowPart) override; //Pointset void RenderWindowPartDeactivated(mitk::IRenderWindowPart* renderWindowPart) override; //Pointset void OnAddPointSetClicked();//Pointset void CreateStreamline(); void RemovefromBundle( bool checked ); void ExtractRandomFibersFromTractogram(); void StartAlgorithm(); void CreatePredictionNode(); void CreateUncertaintySampleNode(); void CreateDistanceSampleNode(); void RemovefromUncertainty( bool checked ); void RemovefromDistance( bool checked ); void RemovefromPrediction( bool checked ); + void RemovefromSelection( bool checked ); void ResampleTractogram(); void RandomPrototypes(); void SFFPrototypes(); void UpdateGui(); ///< update button activity etc. dpending on current datamanager selection protected: void OnCurrentSelectionChanged(QmitkSingleNodeSelectionWidget::NodeList nodes);//Pointset virtual void OnSelectionChanged(berry::IWorkbenchPart::Pointer part, const QList& nodes) override; void OnEndInteraction(); void CreateStreamlineInteractor(); Ui::QmitkInteractiveFiberDissectionViewControls* m_Controls; int m_IterationCounter; ///< used for data node naming int m_RandomExtractionCounter; ///< used for random extracton of different Streamlines int m_activeCycleCounter; + int m_initRandom; // int m_thresh2; std::vector m_SelectedFB; ///< selected fiber bundle nodes // mitk::DataNode::Pointer m_trainbundle; mitk::Image::Pointer m_SelectedImage; mitk::DataNode::Pointer m_SelectedPS; mitk::DataNode::Pointer m_SelectedImageNode; mitk::FiberBundle::Pointer m_positiveBundle; mitk::FiberBundle::Pointer m_newfibersBundle; mitk::FiberBundle::Pointer m_negativeBundle; mitk::FiberBundle::Pointer m_Prediction; mitk::FiberBundle::Pointer m_UncertaintyLabel; mitk::FiberBundle::Pointer m_DistanceLabel; mitk::DataNode::Pointer m_positivBundlesNode; mitk::DataNode::Pointer m_newfibersBundleNode; mitk::DataNode::Pointer m_negativeBundleNode; mitk::DataNode::Pointer m_PredictionNode; mitk::DataNode::Pointer m_UncertaintyLabelNode; mitk::DataNode::Pointer m_DistanceLabelNode; vtkSmartPointer m_positiveFibersData; vtkSmartPointer m_newfibersFibersData; vtkSmartPointer m_picker1; mitk::StreamlineInteractor::Pointer m_StreamlineInteractor; std::shared_ptr< mitk::StreamlineFeatureExtractor > clusterer; std::vector> m_index; }; #endif // _QMITKFIBERTRACKINGVIEW_H_INCLUDED diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionViewControls.ui b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionViewControls.ui index 6acb62f..b4b89d4 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionViewControls.ui +++ b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkInteractiveFiberDissectionViewControls.ui @@ -1,804 +1,811 @@ QmitkInteractiveFiberDissectionViewControls 0 0 417 711 Form QCommandLinkButton:disabled { border: none; } QGroupBox { background-color: transparent; } 9 9 9 9 Qt::Vertical QSizePolicy::Fixed 20 40 Please Select Input Data Input DTI Fiber Bundle: <html><head/><body><p><span style=" color:#ff0000;">mandatory</span></p></body></html> true true true Label Fibers QFrame::NoFrame - 0 + 3 6 0 0 399 455 Generate Prototypes Use Standard Prototypes 10 200 100 Choose own Prototypes Random Extraction Number of Prototypes true Use Standard Prototypes true SFF Extraction Extract from Qt::Vertical 20 40 Resample Bundle 0 0 399 455 Fiber Creation QFrame::StyledPanel QFrame::Raised Selected Streamline Pointset Qt::Horizontal 40 20 Streamline Points Qt::Vertical 20 40 Qt::Horizontal 40 20 Qt::Horizontal 40 20 Add new Streamline Pointset 0 0 0 40 Qt::Vertical 20 40 Create Streamline true 0 0 438 441 Fiber Labelling Dissect/Eraze Fibers by Erasion and Highlighting Qt::Vertical 20 40 true ArrowCursor QFrame::NoFrame QFrame::Raised Label individual Streamlines Streamlines to be labeled QFrame::StyledPanel QFrame::Raised Qt::Horizontal 40 20 Qt::Horizontal 40 20 Qt::Horizontal 40 20 Add Accept: Rightclick + alt 1 2000 10 Random Streamlines Qt::Vertical 20 40 30 30 :/org.mitk.gui.qt.diffusionimaging.fiberprocessing/resources/eraze.png:/org.mitk.gui.qt.diffusionimaging.fiberprocessing/resources/eraze.png 32 32 true Label Streamlines Qt::Horizontal 40 20 Reject: Rightclick + shift Unlabel: Leftclick + shift/alt 0 0 399 455 Active Learning 9 9 341 391 QFrame::StyledPanel QFrame::Raised 1 500 0 10 Qt::Horizontal 40 20 Uncertainty_Label .. Qt::Horizontal 40 20 Add Pred_Label Qt::Horizontal 84 22 Dist: to label Create Prediction Qt::Vertical 20 40 Create Uncertanty Map Add Train Classifier Unc: to label 1 500 10 Qt::Horizontal 40 20 Unc_dis_Label + + + + Selection_Label + + + + + QmitkDataStorageComboBox + QComboBox +
QmitkDataStorageComboBox.h
+
QmitkSingleNodeSelectionWidget QWidget
QmitkSingleNodeSelectionWidget.h
QmitkPointListWidget QWidget
QmitkPointListWidget.h
- - QmitkDataStorageComboBox - QComboBox -
QmitkDataStorageComboBox.h
-