DGtal  0.9.4.1
ITKReader.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
19  * @author Pierre Gueth (\c pierre.gueth@gmail.com )
20  * Laboratoire d'InfoRmatique en Image et Systèmes d'information - LIRIS (CNRS,
21  * UMR 5205), CNRS, France
22  *
23  *
24  * @author Bertrand Kerautret (\c bertrand.kerautret@loria.fr )
25  * LORIA (CNRS, UMR 7503), University of Lorraine, France
26  *
27  * @date 2013/10/28
28  *
29  * Header file for module ITKReader.cpp
30  *
31  * This file is part of the DGtal library.
32  */
33 
34 #include "DGtal/images/ConstImageAdapter.h"
35 #include "DGtal/images/ImageContainerByITKImage.h"
36 #if defined(__GNUG__)
37 #pragma GCC diagnostic push
38 #pragma GCC diagnostic ignored "-Wdeprecated-declarations"
39 #endif
40 #if defined(__clang__)
41 #pragma clang diagnostic push
42 #pragma clang diagnostic ignored "-Wdocumentation"
43 #endif
44 #include <itkImageFileReader.h>
45 #if defined(__clang__)
46 #pragma clang diagnostic pop
47 #endif
48 #if defined(__GNUG__)
49 #endif
50 #pragma GCC diagnostic pop
51 
52 namespace DGtal {
53 
54  template <typename I>
55  template <typename TFunctor>
56  I ITKReader<I>::importITK(
57  const std::string & filename,
58  const TFunctor & aFunctor )
59  {
60  typedef typename Image::Domain Domain;
61  typedef typename Domain::Point Point;
62  typedef itk::ImageIOBase::IOComponentType IOComponentType;
63  BOOST_CONCEPT_ASSERT( (concepts::CUnaryFunctor<TFunctor, ValueOut, Value>));
64  const IOComponentType componentType = getITKComponentType( filename );
65  typedef unsigned char ElementType;
66  switch ( componentType )
67  {
68 
69  default:
70  case itk::ImageIOBase::UNKNOWNCOMPONENTTYPE:
71  std::cerr << "Unknown and unsupported component type!" << std::endl;
72  case itk::ImageIOBase::UCHAR:
73  {
74  typedef ImageContainerByITKImage<Domain, unsigned char> DGtalITKImage;
75  Domain d;
76  DGtalITKImage im( d );
77  return readDGtalImageFromITKtypes<DGtalITKImage>( filename, aFunctor );
78  }
79  case itk::ImageIOBase::CHAR:
80  {
81  typedef ImageContainerByITKImage<Domain, char> DGtalITKImage;
82  Domain d;
83  DGtalITKImage im( d );
84  return readDGtalImageFromITKtypes<DGtalITKImage>( filename, aFunctor );
85  }
86  case itk::ImageIOBase::USHORT:
87  {
88  typedef ImageContainerByITKImage<Domain, unsigned short> DGtalITKImage;
89  Domain d;
90  DGtalITKImage im( d );
91  return readDGtalImageFromITKtypes<DGtalITKImage>( filename, aFunctor );
92  }
93  case itk::ImageIOBase::SHORT:
94  {
95  typedef ImageContainerByITKImage<Domain, short> DGtalITKImage;
96  Domain d;
97  DGtalITKImage im( d );
98  return readDGtalImageFromITKtypes<DGtalITKImage>( filename, aFunctor );
99  }
100  case itk::ImageIOBase::ULONG:
101  {
102  typedef ImageContainerByITKImage<Domain, unsigned long> DGtalITKImage;
103  Domain d;
104  DGtalITKImage im( d );
105  return readDGtalImageFromITKtypes<DGtalITKImage>( filename, aFunctor );
106  }
107  case itk::ImageIOBase::LONG:
108  {
109  typedef ImageContainerByITKImage<Domain, long> DGtalITKImage;
110  Domain d;
111  DGtalITKImage im( d );
112  return readDGtalImageFromITKtypes<DGtalITKImage>( filename, aFunctor );
113  }
114  case itk::ImageIOBase::FLOAT:
115  {
116  typedef ImageContainerByITKImage<Domain, float> DGtalITKImage;
117  Domain d;
118  DGtalITKImage im( d );
119  return readDGtalImageFromITKtypes<DGtalITKImage>( filename, aFunctor );
120  }
121  case itk::ImageIOBase::DOUBLE:
122  {
123  typedef ImageContainerByITKImage<Domain, float> DGtalITKImage;
124  Domain d;
125  DGtalITKImage im( d );
126  return readDGtalImageFromITKtypes<DGtalITKImage>( filename, aFunctor );
127  }
128  }
129  }
130 
131  template <typename I>
132  itk::ImageIOBase::IOComponentType
133  ITKReader<I>::getITKComponentType( const std::string & filename )
134  {
135  typedef itk::ImageIOBase::IOComponentType IOComponentType;
136  IOComponentType componentType;
137  try
138  {
139  itk::ImageIOBase::Pointer imageIO = itk::ImageIOFactory::CreateImageIO(
140  filename.c_str(), itk::ImageIOFactory::ReadMode );
141  imageIO->SetFileName( filename.c_str() );
142  imageIO->ReadImageInformation();
143  componentType = imageIO->GetComponentType();
144  }
145  catch ( itk::ExceptionObject & e )
146  {
147  trace.error() << e;
148  throw IOException();
149  }
150  return componentType;
151  }
152 
153  template <typename I>
154  template <typename TypeDGtalImage, typename TFunctor>
155  typename ITKReader<I>::Image ITKReader<I>::readDGtalImageFromITKtypes(
156  const std::string & filename,
157  const TFunctor & aFunctor )
158  {
159 
160  typedef typename TypeDGtalImage::ITKImagePointer ITKImagePointer;
161  typedef typename Image::Domain Domain;
162 
163  ITKImagePointer itk_image = TypeDGtalImage::ITKImage::New();
164 
165  try
166  {
167  typedef itk::ImageFileReader<typename TypeDGtalImage::ITKImage>
168  ITKImageReader;
169  typename ITKImageReader::Pointer reader = ITKImageReader::New();
170 
171  reader->SetFileName( filename );
172  reader->Update();
173  reader->GetOutput();
174 
175  itk_image = reader->GetOutput();
176  }
177  catch (itk::ExceptionObject &e)
178  {
179  trace.error() << e;
180  throw IOException();
181  }
182 
183  const TypeDGtalImage dgtal_itk_image( itk_image );
184  const Domain& domain = dgtal_itk_image.domain();
185 
186  typedef ConstImageAdapter<TypeDGtalImage, Domain, functors::Identity, Value,
187  TFunctor>
188  AdaptedImage;
189  const functors::Identity identityFunctor{};
190  const AdaptedImage adapted(dgtal_itk_image, domain, identityFunctor, aFunctor);
191 
192  Image image(domain);
193  std::copy(adapted.constRange().begin(), adapted.constRange().end(), image.range().outputIterator());
194  return image;
195  }
196 
197 }//namespace