3dImageViewer.cpp

 
#include <iostream>
 
#include "DGtal/base/Common.h"
#include "DGtal/base/BasicFunctors.h"
#include "DGtal/helpers/StdDefs.h"
#include "DGtal/io/readers/VolReader.h"
#include "DGtal/io/viewers/PolyscopeViewer.h"
#include "DGtal/io/readers/PointListReader.h"
#include "DGtal/io/readers/MeshReader.h"
#include "DGtal/topology/helpers/Surfaces.h"
#include "DGtal/topology/SurfelAdjacency.h"
 
#include "DGtal/io/Color.h"
#include "DGtal/io/colormaps/GradientColorMap.h"
#include "DGtal/io/readers/GenericReader.h"
#ifdef DGTAL_WITH_ITK
#include "DGtal/io/readers/DicomReader.h"
#endif
 
#include "DGtal/images/ImageSelector.h"
 
 
 
// #include "specificClasses/Viewer3DImage.cpp"
 
#include "CLI11.hpp"
 
 
 
using namespace std;
using namespace DGtal;
using namespace Z3i;
 
 
int main( int argc, char** argv )
{
 
  typedef DGtal::ImageContainerBySTLVector<DGtal::Z3i::Domain,  unsigned char > Image3D;
  typedef DGtal::ImageContainerBySTLVector<DGtal::Z2i::Domain,  unsigned char > Image2D;
 
  // parse command line using CLI ----------------------------------------------
  CLI::App app;
  std::string inputFileName;
  std::string inputFileNameSDP;
  std::string inputFileNameMesh;
 
  DGtal::int64_t rescaleInputMin {0};
  DGtal::int64_t rescaleInputMax {255};
  bool grid {false};
  bool intergrid {false};
  bool emptyMode {false};
  bool displayDigitalSurface {false};
  bool thresholdImage {false};
  bool colorizeCC {false};
  int thresholdMin {0};
  int thresholdMax {255};
  std::vector<unsigned int> vectSDPIndex {0,1,2};
  std::vector<unsigned int> colorSDP;
  std::vector<unsigned int> colorMesh;
 
  float sx {1.0};
  float sy {1.0};
  float sz {1.0};
  double ballRadius = {0.0};
  unsigned char transp {255};
 
  
  app.description("Displays volume file as a voxel set by using PolyscopeViewer\n 3dImageViewer  $DGtal/examples/samples/lobster.vol --thresholdImage -m 180");
  
  app.add_option("-i,--input,1", inputFileName, "vol file (.vol, .longvol .p3d, .pgm3d and if DGTAL_WITH_ITK is selected: dicom, dcm, mha, mhd). For longvol, dicom, dcm, mha or mhd formats, the input values are linearly scaled between 0 and 255." )
  ->required()
  ->check(CLI::ExistingFile);
  
  app.add_flag("--grid", grid , "draw slice images using grid mode.");
  app.add_flag("--intergrid", grid , "draw slice images using inter grid mode.");
  app.add_flag("--emptyMode", emptyMode,"remove the default boundingbox display.");
  app.add_flag("--thresholdImage", thresholdImage,"threshold the image to define binary shape");
  app.add_option("--thresholdMin", thresholdMin, "threshold min to define binary shape");
  app.add_option("--thresholdMax", thresholdMax, "threshold maw to define binary shape");
  app.add_option("--displaySDP,s",inputFileNameSDP, "display a set of discrete points (.sdp)" );
  app.add_option("--SDPindex", vectSDPIndex, "specify the sdp index.")
   ->expected(3);
  app.add_option("--SDPball",ballRadius, "use balls to display a set of discrete points (if not set to 0 and used with displaySDP option).");
  
  app.add_option("--displayMesh", inputFileNameMesh, "display a Mesh given in OFF or OFS format.");
  app.add_flag("--displayDigitalSurface",displayDigitalSurface, "display the digital surface instead of display all the set of voxels (used with thresholdImage or displaySDP options)" );
  
  app.add_flag("--colorizeCC", colorizeCC, "colorize each Connected Components of the surface displayed by displayDigitalSurface option.");
  app.add_option("--colorSDP,-c", colorSDP, "set the color  discrete points: r g b a ")
   ->expected(4);
  app.add_option("--colorMesh", colorMesh, "set the color of Mesh (given from displayMesh option) : r g b a ")
   ->expected(4);
 
  app.add_option("--scaleX,-x", sx, "set the scale value in the X direction" );
  app.add_option("--scaleY,-y", sy, "set the scale value in the Y direction" );
  app.add_option("--scaleZ,-z", sy, "set the scale value in the Z direction" );
  app.add_option("--rescaleInputMin",rescaleInputMin, "min value used to rescale the input intensity (to avoid basic cast into 8  bits image)." );
  app.add_option("--rescaleInputMax",rescaleInputMax, "max value used to rescale the input intensity (to avoid basic cast into 8  bits image)." );
  app.add_option("--transparency,-t",transp, "change the default transparency" );
 
  
  app.get_formatter()->column_width(40);
  CLI11_PARSE(app, argc, argv);
  // END parse command line using CLI ----------------------------------------------
 
  
  string extension = inputFileName.substr(inputFileName.find_last_of(".") + 1);
 
  PolyscopeViewer<> viewer;
 
  typedef DGtal::functors::Rescaling<DGtal::int64_t ,unsigned char > RescalFCT;
  Image3D image =  GenericReader< Image3D >::importWithValueFunctor( inputFileName,RescalFCT(rescaleInputMin,
                                                                                             rescaleInputMax,
                                                                                             0,
                                                                                             255));
  Domain domain = image.domain();
 
  trace.info() << "Image loaded: "<<image<< std::endl;
  viewer << image;
  
  // Used to display 3D surface
  Z3i::DigitalSet set3d(domain);
 
  if(thresholdImage){
    viewer.newCubeList("Threshold image");
    viewer.allowReuseList = true;
    for(Domain::ConstIterator it = domain.begin(), itend=domain.end(); it!=itend; ++it){
      unsigned char  val= image( (*it) );
      if(val<=thresholdMax && val >=thresholdMin)
      {
        if(!displayDigitalSurface)
        {
          viewer << WithQuantity(*it, "value", val);
        }
      }else
      {
        set3d.insert(*it);
      }
    }
    viewer.endCurrentGroup();
  }
 
  if(inputFileNameSDP != "" ){
    if(colorSDP.size()==4){
      Color c(colorSDP[0], colorSDP[1], colorSDP[2], colorSDP[3]);
      viewer << c;
    }
 
    vector<Z3i::Point> vectVoxels;
    if(vectSDPIndex.size()==3)
    {
      vectVoxels = PointListReader<Z3i::Point>::getPointsFromFile(inputFileNameSDP, vectSDPIndex);
    }else
    {
      vectVoxels = PointListReader<Z3i::Point>::getPointsFromFile(inputFileNameSDP);
    }
 
    if (ballRadius != 0.0) viewer.drawAsBalls();
    for(unsigned int i=0;i< vectVoxels.size(); i++)
    {
      if(!displayDigitalSurface)
      {
          viewer << vectVoxels.at(i);
      }
      else
      {
        set3d.insert(vectVoxels.at(i));
      }
    }
    viewer.drawAsPaving();
  }
 
  if(inputFileNameMesh != "")
  {
    if(colorMesh.size() != 0)
    {
      Color c(colorMesh[0], colorMesh[1], colorMesh[2], colorMesh[3]);
      viewer.drawColor(c);
    }
    DGtal::Mesh<Z3i::RealPoint> aMesh(colorMesh.size() == 0);
    MeshReader<Z3i::RealPoint>::importOFFFile(inputFileNameMesh, aMesh);
    viewer << aMesh;
  }
 
  if(displayDigitalSurface)
  {
    KSpace K;
    Point low = domain.lowerBound(); low[0]=low[0]-1; low[1]=low[1]-1; low[2]=low[2]-1;
    Point upp = domain.upperBound(); upp[0]=upp[0]+1; upp[1]=upp[1]+1; upp[2]=upp[2]+1;
    K.init(low, upp , true);
    SurfelAdjacency<3> SAdj( true );
    vector<vector<SCell> > vectConnectedSCell;
    trace.info() << "Extracting surface  set ... " ;
    Surfaces<KSpace>::extractAllConnectedSCell(vectConnectedSCell,K, SAdj, set3d, true);
    trace.info()<< " [done] " <<std::endl;
 
    viewer.drawAsSimplified();
    for(unsigned int i= 0; i <vectConnectedSCell.size(); i++)
    {
      for(unsigned int j= 0; j <vectConnectedSCell.at(i).size(); j++)
      {
        const auto& toDraw = vectConnectedSCell.at(i).at(j);
        if(colorizeCC) 
        {
          viewer << WithQuantity(toDraw, "index", i);
        } 
        else if(colorSDP.size() != 0)
        {
          Color c(colorSDP[0], colorSDP[1], colorSDP[2], colorSDP[3]);
          viewer << WithQuantity(toDraw, "color", c);
        } 
        else 
        {
          viewer << toDraw; 
        }
      }
    }
  }
  
  DGtal::Z3i::Point size = image.domain().upperBound() - image.domain().lowerBound();
  DGtal::Z3i::Point center = image.domain().lowerBound()+size/2;
  unsigned int maxDist = std::max(std::max(size[2], size[1]), size[0]);
  viewer.show();
  return 0;
}