DGtal  1.4.beta
DicomReader.ih
1 /**
2  * This program is free software: you can redistribute it and/or modify
3  * it under the terms of the GNU Lesser General Public License as
4  * published by the Free Software Foundation, either version 3 of the
5  * License, or (at your option) any later version.
6  *
7  * This program is distributed in the hope that it will be useful,
8  * but WITHOUT ANY WARRANTY; without even the implied warranty of
9  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
10  * GNU General Public License for more details.
11  *
12  * You should have received a copy of the GNU General Public License
13  * along with this program. If not, see <http://www.gnu.org/licenses/>.
14  *
15  **/
16 
17 /**
18  * @file DicomReader.ih
19  * @author Adrien Krähenbühl (\c adrien.krahenbuhl@loria.fr )
20  * LORIA (CNRS, UMR 7503), Université de Lorraine, France
21  *
22  * @date 2013/10/10
23  *
24  * Implementation of inline methods defined in DicomReader.h
25  *
26  * This file is part of the DGtal library.
27  */
28 
29 ///////////////////////////////////////////////////////////////////////////////
30 // IMPLEMENTATION of inline methods.
31 ///////////////////////////////////////////////////////////////////////////////
32 
33 //////////////////////////////////////////////////////////////////////////////
34 #include <cstdlib>
35 #include <iostream>
36 #include <fstream>
37 #include <sstream>
38 
39 #include "DGtal/io/Color.h"
40 #include "DGtal/images/ConstImageAdapter.h"
41 
42 // Required ITK files to read serie DICOM files
43 // DGtal must be compiled with "-DWITH_ITK=true" option
44 #if defined(__GNUG__)
45 #pragma GCC diagnostic push
46 #pragma GCC diagnostic ignored "-Wdeprecated-declarations"
47 #pragma GCC diagnostic ignored "-Wpedantic"
48 #if __GNUC__ >= 9
49 #pragma GCC diagnostic ignored "-Wdeprecated-copy"
50 #endif
51 #endif
52 #if defined(__clang__)
53 #pragma clang diagnostic push
54 #pragma clang diagnostic ignored "-Wdocumentation"
55 #endif
56 #include <itkImage.h>
57 #include <itkImageSeriesReader.h>
58 #include <itkGDCMImageIO.h>
59 #include <itkGDCMSeriesFileNames.h>
60 #if defined(__clang__)
61 #pragma clang diagnostic pop
62 #endif
63 #if defined(__GNUG__)
64 #pragma GCC diagnostic pop
65 #endif
66 //////////////////////////////////////////////////////////////////////////////
67 
68 
69 
70 ///////////////////////////////////////////////////////////////////////////////
71 // Implementation of inline methods //
72 
73 
74 ///////////////////////////////////////////////////////////////////////////////
75 // Implementation of inline functions and external operators //
76 
77 
78 namespace DGtal {
79 
80 template <typename TImageContainer, typename TFunctor>
81 template <typename Domain, typename PixelType>
82 inline
83 ImageContainerByITKImage<Domain, PixelType>
84 DicomReader<TImageContainer, TFunctor>::
85 importDicomFiles_( const std::vector<std::string> & filenames )
86 {
87  // Definition of image type
88  const unsigned int dimension = Domain::dimension;
89  typedef itk::Image<PixelType, dimension> ItkImage;
90  typedef itk::ImageSeriesReader<ItkImage> ItkReader;
91  typename ItkReader::Pointer reader = ItkReader::New();
92 
93  // Definition of ITK Dicom reader
94  //typedef itk::GDCMImageIO ItkImageIO;
95  //ItkImageIO::Pointer dicomIO = ItkImageIO::New();
96  //reader->SetImageIO( dicomIO );
97 
98  // Series reader
99  reader->SetFileNames( filenames );
100 
101  typedef ImageContainerByITKImage<Domain, PixelType> TypeDGtalImage;
102  typedef typename TypeDGtalImage::ITKImagePointer ITKImagePointer;
103  ITKImagePointer itkImage = nullptr;
104 
105  // Image reading
106  try
107  {
108  reader->Update();
109  }
110  catch( ... )
111  {
112  trace.error() << "DicomReader: can't read " << filenames.size()
113  << " files"<<std::endl;
114  throw IOException();
115  }
116 
117  itkImage = reader->GetOutput();
118 
119  const typename ItkImage::SizeType& inputSize =
120  itkImage->GetLargestPossibleRegion().GetSize();
121  const auto width = inputSize[0];
122  const auto height = inputSize[1];
123  const auto depth = inputSize[2];
124  if ( !height || !width || !depth )
125  {
126  trace.error() << "DicomReader: one dimension is null (w=" << width
127  << ", h=" << height << ", d=" << depth << ")"
128  << std::endl;
129  throw IOException();
130  }
131 
132  const TypeDGtalImage dgtalItkImage( itkImage );
133 
134  return dgtalItkImage;
135 }
136 
137 template <typename TImageContainer, typename TFunctor>
138 template <typename Image, typename Domain, typename OutPixelType,
139  typename PixelType>
140 Image
141 DicomReader<TImageContainer, TFunctor>::Aux<Image, Domain, OutPixelType,
142  PixelType>::
143 importDicomFiles( const std::vector<std::string> & filenames,
144  const TFunctor & aFunctor )
145 {
146  typedef OutPixelType Value;
147 
148  typedef ImageContainerByITKImage<Domain, PixelType> TypeDGtalImage;
149  const TypeDGtalImage dgtalItkImage =
150  importDicomFiles_<Domain, PixelType>( filenames );
151 
152  const Domain& domain = dgtalItkImage.domain();
153 
154  typedef ConstImageAdapter<TypeDGtalImage, Domain, functors::Identity,
155  Value, TFunctor> AdaptedImage;
156  const functors::Identity identityFunctor{};
157  const AdaptedImage adapted( dgtalItkImage, domain, identityFunctor,
158  aFunctor );
159 
160  TImageContainer image( domain );
161  std::copy( adapted.constRange().begin(), adapted.constRange().end(),
162  image.range().outputIterator() );
163  return image;
164 }
165 
166 
167 //specialization
168 template <typename TImageContainer, typename TFunctor>
169 template <typename Domain, typename OutPixelType, typename PixelType>
170 ImageContainerByITKImage<Domain, OutPixelType>
171 DicomReader<TImageContainer, TFunctor>::
172 Aux<ImageContainerByITKImage<Domain, OutPixelType>, Domain, OutPixelType,
173  PixelType>::
174 importDicomFiles( const std::vector<std::string> & filenames,
175  const TFunctor & aFunctor )
176 {
177  typedef OutPixelType Value;
178 
179  typedef ImageContainerByITKImage<Domain, PixelType> TypeDGtalImage;
180  const TypeDGtalImage dgtalItkImage =
181  importDicomFiles_<Domain, PixelType>( filenames );
182 
183  const Domain& domain = dgtalItkImage.domain();
184 
185  typedef ConstImageAdapter<TypeDGtalImage, Domain, functors::Identity, Value,
186  TFunctor> AdaptedImage;
187  const functors::Identity identityFunctor{};
188  const AdaptedImage adapted( dgtalItkImage, domain, identityFunctor, aFunctor );
189 
190  TImageContainer image( domain );
191  std::copy( adapted.constRange().begin(), adapted.constRange().end(),
192  image.range().outputIterator() );
193 
194  //copy ITKImage parameters
195  image.getITKImagePointer()->SetOrigin(
196  dgtalItkImage.getITKImagePointer()->GetOrigin() );
197  image.getITKImagePointer()->SetSpacing(
198  dgtalItkImage.getITKImagePointer()->GetSpacing() );
199  image.getITKImagePointer()->SetDirection(
200  dgtalItkImage.getITKImagePointer()->GetDirection() );
201 
202  return image;
203 }
204 
205 
206 template <typename TImageContainer, typename TFunctor>
207 template <typename PixelType>
208 inline TImageContainer
209 DicomReader<TImageContainer, TFunctor>::
210 importDicomFiles( const std::vector<std::string> & filenames,
211  const TFunctor & aFunctor )
212 {
213  typedef typename TImageContainer::Domain Domain;
214  return Aux<TImageContainer, Domain, Value, PixelType>::
215  importDicomFiles( filenames, aFunctor );
216 }
217 
218 
219 template <typename TImageContainer, typename TFunctor>
220 inline
221 TImageContainer
222 DicomReader<TImageContainer, TFunctor>::
223 importDicom( const std::string & aFilename,
224  const TFunctor & aFunctor )
225 {
226  std::string directory = aFilename.substr( 0, aFilename.find_last_of("/") );
227 
228  typedef itk::GDCMSeriesFileNames ItkNamesGenerator;
229  ItkNamesGenerator::Pointer nameGenerator = ItkNamesGenerator::New();
230  nameGenerator->SetDirectory( directory );
231  const std::vector<std::string> &filenames =
232  nameGenerator->GetInputFileNames();
233 
234  return importDicomFiles<PixelType>( filenames, aFunctor );
235 }
236 
237 template <typename TImageContainer, typename TFunctor>
238 inline
239 TImageContainer
240 DicomReader<TImageContainer, TFunctor>::
241 importDicomSeries( const std::vector<std::string> & filenames,
242  const Functor & aFunctor )
243 {
244  return importDicomFiles<PixelType>( filenames, aFunctor );
245 }
246 
247 }
248 
249 // //
250 ///////////////////////////////////////////////////////////////////////////////
251 
252