DGtalTools  1.2.0
volMask.cpp
1 
32 #include "DGtal/base/Common.h"
33 #include <DGtal/helpers/StdDefs.h>
34 #include <DGtal/images/ImageContainerBySTLVector.h>
35 
36 #include <DGtal/io/readers/GenericReader.h>
37 #include <DGtal/io/writers/GenericWriter.h>
38 
39 #include "CLI11.hpp"
40 
41 #ifdef WITH_ITK
42 #include "DGtal/io/readers/ITKReader.h"
43 #include "DGtal/io/writers/ITKWriter.h"
44 
45 #endif
46 
48 using namespace std;
49 using namespace DGtal;
51 
92 typedef ImageContainerBySTLVector < Z3i::Domain, unsigned char > Image3D;
93 #ifdef WITH_ITK
94 typedef ImageContainerBySTLVector < Z3i::Domain, double > Image3D_D;
95 typedef ImageContainerBySTLVector < Z3i::Domain, int > Image3D_I;
96 #endif
97 
102 template<typename TImage, typename TImageMask>
103 typename TImageMask::Domain
104 subDomainMasked(const TImage &image, const TImageMask &maskImage,
105  typename TImageMask::Value maskValue, unsigned int domainOffset=1){
106  typename TImageMask::Domain res;
107  Z3i::Point minP = image.domain().upperBound();
108  Z3i::Point maxP = image.domain().lowerBound();
109 
110  Z3i::Point::Iterator minIt;
111  Z3i::Point::Iterator maxIt;
112  bool foundMaskedVal = false;
113  for(const auto &p: image.domain())
114  {
115  minIt = minP.begin();
116  maxIt = maxP.begin();
117  if( maskImage.domain().isInside(p) && maskImage(p) == maskValue ) // no noise on mask image
118  {
119  foundMaskedVal=true;
120  for(auto pIt=p.begin(); pIt!=p.end();pIt++ )
121  {
122  if( *pIt < *minIt ){*minIt = *pIt;}
123  if( *pIt > *maxIt ){*maxIt = *pIt;}
124  minIt++;
125  maxIt++;
126  }
127  }
128  }
129  if (!foundMaskedVal)
130  {
131  trace.info() << "No masked value found resulting image will be empty." << std::endl;
132  return image.domain();
133 
134  }
135  // offset to avoid problems on borders
136  Z3i::Point offset(domainOffset,domainOffset,domainOffset);
137  minP -= offset;
138  maxP += offset;
139  trace.info() << "sub-domain:" << minP << " " << maxP << std::endl;
140  return typename TImageMask::Domain(minP, maxP);
141 }
142 
143 
144 template<typename TImage, typename TImageMask>
145 void
146 applyMask(const TImage &inputImage,TImage &outputImage,
147  const TImageMask &maskImage, typename TImageMask::Value maskValue)
148 {
149  for (const auto &p: outputImage.domain())
150  {
151  if (inputImage.domain().isInside(p) && maskImage(p) == maskValue)
152  {
153  outputImage.setValue(p, inputImage(p) );
154  }
155  }
156 
157 }
158 
159 template<typename TImage, typename TImageMask>
160 void
161 processImage(const TImage &inputImage, const TImageMask &maskImage,
162  typename TImageMask::Value maskValue, std::string outputFileName, unsigned int offsetBorder ){
163  // First step getting the bounding box of the domain:
164  auto subDm = subDomainMasked(inputImage, maskImage, maskValue, offsetBorder);
165  TImage outputImage( subDm );
166  // Second step: masking source image
167  applyMask(inputImage, outputImage, maskImage, maskValue);
168  trace.info() << "writing output image...";
169  GenericWriter<TImage>::exportFile(outputFileName, outputImage);
170 }
171 
172 
173 int main( int argc, char** argv )
174 {
175 
176  // parse command line using CLI ----------------------------------------------
177  CLI::App app;
178  std::string inputFileName;
179  std::string outputFileName {"result.vol"};
180  std::string inputType;
181  std::string maskFileName;
182  unsigned int offsetBorder {0};
183  int maskValue {1};
184 
185  app.description("Outputs a new image from two input images, one representing the data, one representing the selection mask. The size of output image is the size of the bounding box of selected values, plus the chosen border offset. \n Typical use example:\n \t volMask ${DGtal}/examples/samples/lobster.vol lobsMasked.vol -a ${DGtal}/examples/samples/lobster.vol -m 100 \n");
186 
187 #ifdef WITH_ITK
188  app.add_option("-i,--input,1", inputFileName, "an input 3D image vol (or ITK: .nii, mha, ... ) file." )
189  ->required()
190  ->check(CLI::ExistingFile);
191  app.add_option("--inputType,-t",inputType, "to specify the input image type (int or double).");
192 
193  #else
194  app.add_option("-i,--input,1", inputFileName, "an input vol file." )
195  ->required()
196  ->check(CLI::ExistingFile);
197 #endif
198 
199  app.add_option("--mask,-a",maskFileName, "the mask image that represents the elements that are copied as output in the resulting image (by default set to 1 you can change this value by using --maskValue). ")
200  ->check(CLI::ExistingFile);
201  app.add_option("-o,--output,2", outputFileName, "the output masked image.", true );
202  app.add_option("--offsetBorder,-f", offsetBorder, "add a border offset to the bounding box of the masked value domain.", true);
203  app.add_option("--maskValue,-m", maskValue, "the masking value.", true);
204  app.get_formatter()->column_width(40);
205  CLI11_PARSE(app, argc, argv);
206  // END parse command line using CLI ----------------------------------------------
207 
208 
209 
210  trace.info() << "Reading mask image...";
211  Image3D maskImage = DGtal::GenericReader<Image3D>::import(maskFileName);
212  trace.info() << "[done]"<< std::endl;
213  trace.info() << "Reading input image...";
214 #ifdef WITH_ITK
215  if (inputType=="double")
216  {
217  Image3D_D inputImage = DGtal::GenericReader<Image3D_D>::import(inputFileName);
218  trace.info() << "[done]"<< std::endl;
219  processImage(inputImage, maskImage, maskValue, outputFileName, offsetBorder);
220  }
221  else if (inputType=="int")
222  {
223  Image3D_I inputImage = DGtal::GenericReader<Image3D_I>::import(inputFileName);
224  trace.info() << "[done]"<< std::endl;
225  processImage(inputImage, maskImage, maskValue, outputFileName, offsetBorder);
226  }
227  else
228  {
229  Image3D inputImage = DGtal::GenericReader<Image3D>::import(inputFileName);
230  trace.info() << "[done]"<< std::endl;
231  processImage(inputImage, maskImage, maskValue, outputFileName, offsetBorder);
232  }
233 #else
234  Image3D inputImage = DGtal::GenericReader<Image3D>::import(inputFileName);
235  processImage(inputImage, maskImage, maskValue, outputFileName, offsetBorder);
236 #endif
237 
238  trace.info() << "[Done]" << std::endl;
239  return 0;
240 }
241 
Definition: ATu0v1.h:57