NorMIT-nav  16.5
An IGT application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
cxUSFrameData.cpp
Go to the documentation of this file.
1 /*=========================================================================
2 This file is part of CustusX, an Image Guided Therapy Application.
3 
4 Copyright (c) 2008-2014, SINTEF Department of Medical Technology
5 All rights reserved.
6 
7 Redistribution and use in source and binary forms, with or without
8 modification, are permitted provided that the following conditions are met:
9 
10 1. Redistributions of source code must retain the above copyright notice,
11  this list of conditions and the following disclaimer.
12 
13 2. Redistributions in binary form must reproduce the above copyright notice,
14  this list of conditions and the following disclaimer in the documentation
15  and/or other materials provided with the distribution.
16 
17 3. Neither the name of the copyright holder nor the names of its contributors
18  may be used to endorse or promote products derived from this software
19  without specific prior written permission.
20 
21 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
22 AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
23 IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
24 DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
25 FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
26 DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
27 SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28 CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
29 OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
30 OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31 =========================================================================*/
32 
33 
34 #include "cxUSFrameData.h"
35 
36 #include <vtkImageData.h>
37 #include <vtkImageLuminance.h>
38 #include <vtkImageClip.h>
39 #include <vtkImageAppend.h>
40 #include <vtkMetaImageWriter.h>
41 #include <vtkImageImport.h>
42 #include "cxTypeConversions.h"
43 #include "cxDataReaderWriter.h"
44 #include <QFileInfo>
45 #include "cxTimeKeeper.h"
46 #include "cxImageDataContainer.h"
47 #include "cxVolumeHelpers.h"
48 #include "cxLogger.h"
49 
50 
51 typedef vtkSmartPointer<vtkImageAppend> vtkImageAppendPtr;
52 
53 namespace cx
54 {
55 
56 ProcessedUSInputData::ProcessedUSInputData(std::vector<vtkImageDataPtr> frames, std::vector<TimedPosition> pos, vtkImageDataPtr mask, QString path, QString uid) :
57  mProcessedImage(frames),
58  mFrames(pos),
59  mMask(mask),
60  mPath(path),
61  mUid(uid)
62 {
63  this->validate();
64 }
65 
67 {
68  std::vector<TimedPosition> frameInfo = this->getFrames();
69  Eigen::Array3i inputDims = this->getDimensions();
70 
71  if (inputDims[2] != static_cast<int>(frameInfo.size()))
72  {
73  QString msg("Mismatch in US input data: inputDims[2] != frameInfo.size() : %1 != %2");
74  reportWarning(msg.arg(inputDims[2]).arg(frameInfo.size()));
75  return false;
76  }
77 
78  if (inputDims[2]==0)
79  {
80  reportWarning("Empty US input data");
81  return false;
82  }
83 
84  return true;
85 }
86 
87 unsigned char* ProcessedUSInputData::getFrame(unsigned int index) const
88 {
89  CX_ASSERT(index < mProcessedImage.size());
90 
91  // Raw data pointer
92  unsigned char *inputPointer = static_cast<unsigned char*> (mProcessedImage[index]->GetScalarPointer());
93  return inputPointer;
94 }
95 
97 {
98  Eigen::Array3i retval;
99  retval[0] = mProcessedImage[0]->GetDimensions()[0];
100  retval[1] = mProcessedImage[0]->GetDimensions()[1];
101  retval[2] = mProcessedImage.size();
102  return retval;
103 }
104 
106 {
107  Vector3D retval = Vector3D(mProcessedImage[0]->GetSpacing());
108  retval[2] = retval[0]; // set z-spacing to arbitrary value.
109  return retval;
110 }
111 
112 std::vector<TimedPosition> ProcessedUSInputData::getFrames() const
113 {
114  return mFrames;
115 }
116 
118 {
119  return mMask;
120 }
121 
123 {
124  return mPath;
125 }
126 
128 {
129  return mUid;
130 }
131 
132 
136 
138 {
139  USFrameDataPtr retval(new USFrameData());
140 
141  retval->mName = QFileInfo(inputFrameData->getName()).completeBaseName();
142  retval->mImageContainer.reset(new cx::SplitFramesContainer(inputFrameData->getBaseVtkImageData()));
143  retval->resetRemovedFrames();
144 
145  return retval;
146 }
147 
154 USFrameDataPtr USFrameData::create(QString inputFilename)
155 {
156  QFileInfo info(inputFilename);
157 
158  TimeKeeper timer;
159  QString mhdSingleFile = info.absolutePath()+"/"+info.completeBaseName()+".mhd";
160 
161  if (QFileInfo(mhdSingleFile).exists())
162  {
163  vtkImageDataPtr image = MetaImageReader().loadVtkImageData(mhdSingleFile);
164  // load from single file
165  USFrameDataPtr retval = USFrameData::create(ImagePtr(new Image(mhdSingleFile, image)));
166  retval->mName = QFileInfo(mhdSingleFile).completeBaseName();
167  timer.printElapsedms(QString("Loading single %1").arg(inputFilename));
168  return retval;
169  }
170  else
171  {
172  USFrameDataPtr retval(new USFrameData());
173  retval->mName = QFileInfo(inputFilename).completeBaseName();
174  retval->mImageContainer.reset(new cx::CachedImageDataContainer(inputFilename, -1));
175  retval->resetRemovedFrames();
176  return retval;
177  }
178 }
179 
181 {
182  if (!images)
183  return USFrameDataPtr();
184  USFrameDataPtr retval(new USFrameData());
185  retval->mName = name;
186  retval->mImageContainer = images;
187  retval->resetRemovedFrames();
188 
189  return retval;
190 }
191 
192 USFrameDataPtr USFrameData::create(QString name, std::vector<vtkImageDataPtr> frames)
193 {
194  cx::ImageDataContainerPtr container;
195  container.reset(new cx::FramesDataContainer(frames));
196  return create(name, container);
197 }
198 
199 
201  mCropbox(0,0,0,0,0,0), mPurgeInput(true)
202 {
203 }
204 
206 {
207  if (mImageContainer->empty())
208  return;
209 
210  mReducedToFull.clear();
211  for (unsigned i=0; i<mImageContainer->size(); ++i)
212  mReducedToFull.push_back(i);
213 }
214 
216 {
217  mImageContainer.reset();
218 }
219 
223 void USFrameData::removeFrame(unsigned int index)
224 {
225  CX_ASSERT(index < mReducedToFull.size());
226  mReducedToFull.erase(mReducedToFull.begin()+index);
227 }
228 
229 Eigen::Array3i USFrameData::getDimensions() const
230 {
231  Eigen::Array3i retval(mImageContainer->get(0)->GetDimensions());
232 
233  if (mCropbox.range()[0]!=0)
234  {
235  // WARNING: cannot use dimensions from the cropImage function - it uses extent,
236  // which may be larger. We need the real size, given by wholeextent/updateextent.
237  vtkImageDataPtr sample = this->cropImageExtent(mImageContainer->get(0), mCropbox);
238 // Eigen::Array<int, 6, 1> extent(sample->GetWholeExtent());
239  Eigen::Array<int, 6, 1> extent(sample->GetExtent());
240  retval[0] = extent[1]-extent[0]+1;
241  retval[1] = extent[3]-extent[2]+1;
242 // std::cout << "extent dim: " << retval[0] << " " << retval[1] << std::endl;
243  }
244 
245  retval[2] = mReducedToFull.size();
246  return retval;
247 }
248 
249 
251 {
252  Vector3D retval(1,1,1);
253  if (!mImageContainer->empty())
254  retval = Vector3D(mImageContainer->get(0)->GetSpacing());
255  retval[2] = retval[0]; // set z-spacing to arbitrary value.
256  return retval;
257 }
258 
260 {
261  // ensure clip never happens in z dir.
262  cropbox[4] = -100000;
263  cropbox[5] = 100000;
264  mCropbox = cropbox;
265 }
266 
274 {
275  vtkImageClipPtr clip = vtkImageClipPtr::New();
276  clip->SetInputData(input);
277  clip->SetOutputWholeExtent(cropbox.begin());
278 // clip->ClipDataOn(); // this line sets wholeextent==extent, but uses lots of time (~6ms/frame)
279  clip->Update();
280  vtkImageDataPtr rawResult = clip->GetOutput();
281 // rawResult->Update();
282  rawResult->Crop(cropbox.begin()); // moved in from toGrayscaleAndEffectuateCropping when VTK5->6
283  return rawResult;
284 }
285 
291 {
292  vtkImageDataPtr grayScaleData;
293 
294  if (input->GetNumberOfScalarComponents() == 1) // already gray
295  {
296  // crop (slow)
297  grayScaleData = input;
298 // outData->Crop();
299  }
300  else
301  {
302  // convert and crop as side effect (optimization)
303  grayScaleData = convertImageDataToGrayScale(input);
304  }
305 
306  vtkImageDataPtr outData = this->convertTo8bit(grayScaleData);
307 
308  vtkImageDataPtr copy = vtkImageDataPtr::New();
309  copy->DeepCopy(outData);
310  return copy;
311 }
312 
313 vtkImageDataPtr USFrameData::convertTo8bit(vtkImageDataPtr input) const
314 {
315  vtkImageDataPtr retval = input;
316  if (input->GetScalarSize() > 1)
317  {
318  ImagePtr tempImage = cx::ImagePtr(new cx::Image("tempImage", input, "tempImage"));
319  tempImage->resetTransferFunctions();
320  retval = tempImage->get8bitGrayScaleVtkImageData();
321  }
322  return retval;
323 }
324 
325 // for testing
326 //void printStuff(QString text, vtkImageDataPtr data)
327 //{
328 // std::cout << text << std::endl;
329 // Eigen::Array<int, 6, 1> wholeExtent(data->GetWholeExtent());
330 // Eigen::Array<int, 6, 1> extent(data->GetExtent());
331 // Eigen::Array<int, 6, 1> updateExtent(data->GetUpdateExtent());
332 // std::cout << "\t whole ext " << wholeExtent << std::endl;
333 // std::cout << "\tupdate ext " << updateExtent << std::endl;
334 // std::cout << "\t ext " << extent << std::endl;
336 
337 // Eigen::Array<vtkIdType,3,1> hinc;
338 // data->GetContinuousIncrements(data->GetWholeExtent(), hinc[0], hinc[1], hinc[2]);
339 // std::cout << "\t whole inc " << hinc << std::endl;
340 
341 // Eigen::Array<vtkIdType,3,1> uinc;
342 // data->GetContinuousIncrements(data->GetUpdateExtent(), uinc[0], uinc[1], uinc[2]);
343 // std::cout << "\t uinc " << uinc << std::endl;
344 
345 // Eigen::Array<vtkIdType,3,1> inc;
346 // data->GetContinuousIncrements(data->GetExtent(), inc[0], inc[1], inc[2]);
347 // std::cout << "\t inc " << inc << std::endl;
348 
349 // Eigen::Array3i dim(data->GetDimensions());
350 // std::cout << "\t dim " << dim << std::endl;
351 
352 // Eigen::Array3i max(wholeExtent[1] - wholeExtent[0], wholeExtent[3] - wholeExtent[2], wholeExtent[5] - wholeExtent[4]);
353 // std::cout << "\t max " << max << std::endl;
354 //}
355 
357 {
358  // Some of the code here is borrowed from the vtk examples:
359  // http://public.kitware.com/cgi-bin/viewcvs.cgi/*checkout*/Examples/Build/vtkMy/Imaging/vtkImageFoo.cxx?root=VTK&content-type=text/plain
360 
361  if (inData->GetNumberOfScalarComponents() != 3)
362  {
363  if(frameNum == 0) //Only report warning once
364  reportWarning("Angio requested for grayscale ultrasound");
365  return grayFrame;
366  }
367 
368  vtkImageDataPtr outData = vtkImageDataPtr::New();
369  outData->DeepCopy(grayFrame);
370 // outData->Update(); // updates whole extent.
371 
372 // printStuff("Clipped color in", inData);
373 // printStuff("Grayscale in", outData);
374 
375 // int* inExt = inData->GetWholeExtent();
376  int* outExt = outData->GetExtent();
377 
378  // Remember that the input might (and do due to vtkImageClip) contain leaps.
379  // This means that the wholeextent might be larger than the extent, thus
380  // we must use a startoffset and leaps between lines.
381 
382  unsigned char *inPtr = static_cast<unsigned char*> (inData->GetScalarPointerForExtent(inData->GetExtent()));
383  unsigned char *outPtr = static_cast<unsigned char*> (outData->GetScalarPointerForExtent(outData->GetExtent()));
384 
385  int maxX = outExt[1] - outExt[0];
386  int maxY = outExt[3] - outExt[2];
387  int maxZ = outExt[5] - outExt[4];
388 
389  Eigen::Array<vtkIdType,3,1> inInc;
390  inData->GetContinuousIncrements(inData->GetExtent(), inInc[0], inInc[1], inInc[2]);
391  CX_ASSERT(inInc[0]==0);
392  // we assume (wholeextent == extent) for the outData in the algorithm below. assert here.
393  Eigen::Array<vtkIdType,3,1> outInc;
394  outData->GetContinuousIncrements(outData->GetExtent(), outInc[0], outInc[1], outInc[2]);
395  CX_ASSERT(outInc[0]==0);
396  CX_ASSERT(outInc[1]==0);
397  CX_ASSERT(outInc[2]==0);
398 
399  for (int z=0; z<=maxZ; z++)
400  {
401  for (int y=0; y<=maxY; y++)
402  {
403  for (int x=0; x<=maxX; x++)
404  {
405  //Look at 3 scalar components at the same time (RGB),
406  //if color is grayscale or close to grayscale: set to zero.
407 
408  if (((*inPtr) == (*(inPtr + 1))) && ((*inPtr) == (*(inPtr + 2))))
409  {
410  (*outPtr) = 0;
411  }
412  else
413  {
414  // strong check: look for near-gray values and remove them.
415  double r = inPtr[0];
416  double g = inPtr[1];
417  double b = inPtr[2];
418  int metric = (fabs(r-g) + fabs(r-b) + fabs(g-b)) / 3; // average absolute diff must be less than or equal to this
419 // std::cout << QString(" %1,%2,%3 \t %4, %5 -- %6").arg(int(inPtr[0])).arg(int(inPtr[1])).arg(int(inPtr[2])).arg(idxR).arg(idxY).arg(metric) << std::endl;
420  if (metric <= 3)
421  (*outPtr) = 0;
422 
423  }
424  //Assume the outVolume is treated with the luminance filter first
425  outPtr++;
426  inPtr += 3;
427  }
428  inPtr += inInc[1];
429  }
430  inPtr += inInc[2];
431  }
432 
433  return outData;
434 }
435 
437 {
438  mPurgeInput = value;
439 }
440 
441 std::vector<std::vector<vtkImageDataPtr> > USFrameData::initializeFrames(std::vector<bool> angio)
442 {
443 
444  std::vector<std::vector<vtkImageDataPtr> > raw(angio.size());
445 
446  for (unsigned i=0; i<raw.size(); ++i)
447  {
448  raw[i].resize(mReducedToFull.size());
449  }
450 
451  // apply cropping and angio
452  for (unsigned i=0; i<mReducedToFull.size(); ++i)
453  {
456 
457  if (mCropbox.range()[0]!=0)
458  current = this->cropImageExtent(current, mCropbox);
459 
460  // optimization: grayFrame is used in both calculations: compute once
461  vtkImageDataPtr grayFrame = this->to8bitGrayscaleAndEffectuateCropping(current);
462 
463  for (unsigned j=0; j<angio.size(); ++j)
464  {
465 
466  if (angio[j])
467  {
468  vtkImageDataPtr angioFrame = this->useAngio(current, grayFrame, i);
469  raw[j][i] = angioFrame;
470  }
471  else
472  {
473  raw[j][i] = grayFrame;
474  }
475  }
476 
477  if (mPurgeInput)
478  mImageContainer->purge(mReducedToFull[i]);
479  }
480 
481  if (mPurgeInput)
482  mImageContainer->purgeAll();
483 
484  return raw;
485 }
486 
488 {
489  mImageContainer->purgeAll();
490 }
491 
493 {
494  USFrameDataPtr retval(new USFrameData(*this));
495  this->resetRemovedFrames();
496  return retval;
497 }
498 
499 QString USFrameData::getName() const
500 {
501  return QFileInfo(mName).completeBaseName();
502 }
503 
505 {
506  return mImageContainer->size();
507 }
508 
510 {
511  TimeKeeper timer;
512  vtkImageDataPtr source = mImageContainer->get(index);
513 
514 // use this test code to display angio output:
515 // vtkImageDataPtr current = mImageContainer->get(index);
516 // current = this->cropImage(current, IntBoundingBox3D(157, 747, 68, 680, 0, 0));
517 // vtkImageDataPtr grayFrame = this->toGrayscale(current);
518 // static vtkImageDataPtr source;
519 // source = this->useAngio(current, grayFrame);
520 
521 // source->Update(); // VTK5->6
522 
523  import->SetImportVoidPointer(source->GetScalarPointer());
524  import->SetDataScalarType(source->GetScalarType());
525  import->SetDataSpacing(source->GetSpacing());
526  import->SetNumberOfScalarComponents(source->GetNumberOfScalarComponents());
527  import->SetWholeExtent(source->GetExtent());
528  import->SetDataExtentToWholeExtent();
529 }
530 
532 {
533  int numberOfFrames = mReducedToFull.size();
534 
535  if (numberOfFrames > 0)
536  if(mImageContainer->get(0)->GetDataDimension() == 3)
537  return true;
538  return false;
539 }
540 
542 {
543  if (mImageContainer->get(0)->GetNumberOfScalarComponents() == 1)
544  return true;
545  return false;
546 }
547 
548 }//namespace cx
vtkImageDataPtr cropImageExtent(vtkImageDataPtr input, IntBoundingBox3D cropbox) const
Use only US angio data as input. Removes grayscale from the US data and converts the remaining color ...
bool is8bit() const
unsigned char * getFrame(unsigned int index) const
vtkSmartPointer< vtkImageAppend > vtkImageAppendPtr
std::vector< std::vector< vtkImageDataPtr > > initializeFrames(std::vector< bool > angio)
cx::ImageDataContainerPtr mImageContainer
#define CX_ASSERT(statement)
Definition: cxLogger.h:131
void setPurgeInputDataAfterInitialize(bool value)
boost::shared_ptr< class Image > ImagePtr
Definition: cxDicomWidget.h:48
boost::shared_ptr< class USFrameData > USFrameDataPtr
void printElapsedms(QString text="Elapsed") const
std::vector< int > mReducedToFull
map from indexes in the reduced volume to the full (original) volume.
virtual vtkImageDataPtr loadVtkImageData(QString filename)
Eigen::Array3i getDimensions() const
IntBoundingBox3D mCropbox
unsigned getNumImages()
QString getName() const
std::vector< TimedPosition > getFrames() const
virtual USFrameDataPtr copy()
Eigen::Vector3i range() const
void reportWarning(QString msg)
Definition: cxLogger.cpp:91
vtkImageDataPtr getMask()
Eigen::Array3i getDimensions() const
A volumetric data set.
Definition: cxImage.h:66
vtkImageDataPtr useAngio(vtkImageDataPtr inData, vtkImageDataPtr grayFrame, int frameNum) const
void fillImageImport(vtkImageImportPtr import, int index)
fill import with a single frame
Representation of an integer bounding box in 3D. The data are stored as {xmin,xmax,ymin,ymax,zmin,zmax}, in order to simplify communication with vtk.
vtkSmartPointer< class vtkImageClip > vtkImageClipPtr
Eigen::Vector3d Vector3D
Vector3D is a representation of a point or vector in 3D.
Definition: cxVector3D.h:63
vtkImageDataPtr convertImageDataToGrayScale(vtkImageDataPtr image)
vtkSmartPointer< class vtkImageImport > vtkImageImportPtr
Vector3D getSpacing() const
void setCropBox(IntBoundingBox3D mCropbox)
Vector3D getSpacing() const
vtkSmartPointer< class vtkImageData > vtkImageDataPtr
Reader for metaheader .mhd files.
vtkImageDataPtr to8bitGrayscaleAndEffectuateCropping(vtkImageDataPtr input) const
ProcessedUSInputData(std::vector< vtkImageDataPtr > frames, std::vector< TimedPosition > pos, vtkImageDataPtr mask, QString path, QString uid)
static USFrameDataPtr create(ImagePtr inputFrameData)
void removeFrame(unsigned int index)
boost::shared_ptr< class ImageDataContainer > ImageDataContainerPtr
Definition: cxUSFrameData.h:45