DGtalTools  0.9.4
vol2heightfield.cpp
1 
29 #include <iostream>
31 #include <fstream>
32 #include "DGtal/base/Common.h"
33 #include "DGtal/helpers/StdDefs.h"
34 #include "DGtal/images/ImageContainerBySTLVector.h"
35 #include "DGtal/io/writers/GenericWriter.h"
36 #include "DGtal/io/readers/GenericReader.h"
37 #include "DGtal/images/ConstImageAdapter.h"
38 #include "DGtal/kernel/BasicPointFunctors.h"
39 
40 
41 #include <boost/program_options/options_description.hpp>
42 #include <boost/program_options/parsers.hpp>
43 #include <boost/program_options/variables_map.hpp>
44 
45 using namespace std;
46 using namespace DGtal;
47 
48 
50 namespace po = boost::program_options;
51 
52 
113 int main( int argc, char** argv )
114 {
118  Image3D::Value, DGtal::functors::Identity > ImageAdapterExtractor;
119 
120  // parse command line ----------------------------------------------
121  po::options_description general_opt("Allowed options are: ");
122  general_opt.add_options()
123  ("help,h", "display this message")
124  ("input,i", po::value<std::string>(), "vol file (.vol, .longvol .p3d, .pgm3d and if 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." )
125  ("output,o", po::value<std::string>(), "sequence of discrete point file (.sdp) ")
126  ("thresholdMin,m", po::value<int>()->default_value(128), "min threshold (default 128)" )
127  ("thresholdMax,M", po::value<int>()->default_value(255), "max threshold (default 255)" )
128  ("nx", po::value<double>()->default_value(0), "set the x component of the projection direction." )
129  ("ny", po::value<double>()->default_value(0), "set the y component of the projection direction." )
130  ("nz", po::value<double>()->default_value(1), "set the z component of the projection direction." )
131  ("centerX,x", po::value<unsigned int>()->default_value(0), "choose x center of the projected image." )
132  ("centerY,y", po::value<unsigned int>()->default_value(0), "choose y center of the projected image." )
133  ("centerZ,z", po::value<unsigned int>()->default_value(1), "choose z center of the projected image." )
134  ("width", po::value<unsigned int>()->default_value(100), "set the width of the resulting height Field image." )
135  ("height", po::value<unsigned int>()->default_value(100), "set the height of the resulting height Field image." )
136  ("heightFieldMaxScan", po::value<unsigned int>()->default_value(255), "set the maximal scan deep." )
137  ("setBackgroundLastDepth", "change the default background (black with the last filled intensity).")
138  ("rescaleInputMin", po::value<DGtal::int64_t>()->default_value(0), "min value used to rescale the input intensity (to avoid basic cast into 8 bits image).")
139  ("rescaleInputMax", po::value<DGtal::int64_t>()->default_value(255), "max value used to rescale the input intensity (to avoid basic cast into 8 bits image).");
140 
141 
142 
143  bool parseOK=true;
144  po::variables_map vm;
145  try{
146  po::store(po::parse_command_line(argc, argv, general_opt), vm);
147  }catch(const std::exception& ex){
148  parseOK=false;
149  trace.info()<< "Error checking program options: "<< ex.what()<< endl;
150  }
151  po::notify(vm);
152  if( !parseOK || vm.count("help")||argc<=1)
153  {
154  std::cout << "Usage: " << argv[0] << " [input] [output]\n"
155  << "Convert volumetric file into a projected 2D image given from a normal direction N and from a starting point P. The 3D volume is scanned in this normal direction N starting from P with a step 1. If the intensity of the 3d point is inside the given thresholds its 2D gray values are set to the current scan number."
156  << general_opt << "\n";
157  std::cout << "Example:\n"
158  << "vol2heightfield -i ${DGtal}/examples/samples/lobster.vol -m 60 -M 500 --nx 0 --ny 0.7 --nz -1 -x 150 -y 0 -z 150 --width 300 --height 300 --heightFieldMaxScan 350 -o resultingHeightMap.pgm \n";
159  return 0;
160  }
161 
162  if(! vm.count("input") ||! vm.count("output"))
163  {
164  trace.error() << " Input and output filename are needed to be defined" << endl;
165  return 0;
166  }
167 
168  string inputFilename = vm["input"].as<std::string>();
169  string outputFilename = vm["output"].as<std::string>();
170  DGtal::int64_t rescaleInputMin = vm["rescaleInputMin"].as<DGtal::int64_t>();
171  DGtal::int64_t rescaleInputMax = vm["rescaleInputMax"].as<DGtal::int64_t>();
172 
173  trace.info() << "Reading input file " << inputFilename ;
174 
176  Image3D inputImage = GenericReader< Image3D >::importWithValueFunctor( inputFilename,RescalFCT(rescaleInputMin,
177  rescaleInputMax,
178  0, 255) );
179 
180 
181  trace.info() << " [done] " << std::endl ;
182 
183  std::ofstream outStream;
184  outStream.open(outputFilename.c_str());
185  int minTh = vm["thresholdMin"].as<int>();
186  int maxTh = vm["thresholdMax"].as<int>();
187 
188  trace.info() << "Processing image to output file " << outputFilename ;
189 
190  unsigned int widthImageScan = vm["height"].as<unsigned int>();
191  unsigned int heightImageScan = vm["width"].as<unsigned int>();
192  unsigned int maxScan = vm["heightFieldMaxScan"].as<unsigned int>();
193  if(maxScan > std::numeric_limits<Image2D::Value>::max()){
194  maxScan = std::numeric_limits<Image2D::Value>::max();
195  trace.warning()<< "value --setBackgroundLastDepth outside mox value of image. Set to max value:" << maxScan << std::endl;
196  }
197 
198  unsigned int centerX = vm["centerX"].as<unsigned int>();
199  unsigned int centerY = vm["centerY"].as<unsigned int>();
200  unsigned int centerZ = vm["centerZ"].as<unsigned int>();
201 
202  double nx = vm["nx"].as<double>();
203  double ny = vm["ny"].as<double>();
204  double nz = vm["nz"].as<double>();
205 
206 
207  Image2D::Domain aDomain2D(DGtal::Z2i::Point(0,0),
208  DGtal::Z2i::Point(widthImageScan, heightImageScan));
209  Z3i::Point ptCenter (centerX, centerY, centerZ);
210  Z3i::RealPoint normalDir (nx, ny, nz);
211  Image2D resultingImage(aDomain2D);
212 
213  for(Image2D::Domain::ConstIterator it = resultingImage.domain().begin();
214  it != resultingImage.domain().end(); it++){
215  resultingImage.setValue(*it, 0);
216  }
218 
219  unsigned int maxDepthFound = 0;
220  for(unsigned int k=0; k < maxScan; k++){
222  ptCenter+normalDir*k,
223  normalDir,
224  widthImageScan);
225  ImageAdapterExtractor extractedImage(inputImage, aDomain2D, embedder, idV);
226  for(Image2D::Domain::ConstIterator it = extractedImage.domain().begin();
227  it != extractedImage.domain().end(); it++){
228  if( resultingImage(*it)== 0 && extractedImage(*it) < maxTh &&
229  extractedImage(*it) > minTh){
230  maxDepthFound = k;
231  resultingImage.setValue(*it, maxScan-k);
232  }
233  }
234  }
235  if (vm.count("setBackgroundLastDepth")){
236  for(Image2D::Domain::ConstIterator it = resultingImage.domain().begin();
237  it != resultingImage.domain().end(); it++){
238  if( resultingImage(*it)== 0 ){
239  resultingImage.setValue(*it, maxScan-maxDepthFound);
240  }
241  }
242  }
243 
244  resultingImage >> outputFilename;
245 
246  trace.info() << " [done] " << std::endl ;
247  return 0;
248 }
249 
250 
251 
252 
STL namespace.
Trace trace(traceWriterTerm)
std::ostream & warning()
std::ostream & info()
const Domain & domain() const
std::vector< Value >::const_iterator ConstIterator
std::ostream & error()
boost::int64_t int64_t