DGtal  1.2.0
VolReader.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 VolReader.ih
19  * @author David Coeurjolly (\c david.coeurjolly@liris.cnrs.fr )
20  * Laboratoire d'InfoRmatique en Image et Systèmes d'information - LIRIS (CNRS, UMR 5205), CNRS, France
21  *
22  * @date 2010/07/25
23  *
24  * Implementation of inline methods defined in VolReader.h
25  *
26  * This file is part of the DGtal library.
27  */
28 
29 
30 //////////////////////////////////////////////////////////////////////////////
31 #include <cstdlib>
32 #include <boost/iostreams/filtering_streambuf.hpp>
33 #include <boost/iostreams/copy.hpp>
34 #include <boost/iostreams/filter/zlib.hpp>
35 //////////////////////////////////////////////////////////////////////////////
36 
37 
38 ///////////////////////////////////////////////////////////////////////////////
39 // Interface - public :
40 
41 
42 #define MAX_HEADERNUMLINES 64
43 
44 
45 
46 template <typename T, typename TFunctor>
47 inline
48 T
49 DGtal::VolReader<T, TFunctor>::importVol( const std::string & filename,
50  const Functor & aFunctor)
51 {
52  FILE * fin;
53  DGtal::IOException dgtalexception;
54 
55 
56  typename T::Point firstPoint( 0, 0, 0 );
57  typename T::Point lastPoint( 0, 0, 0 );
58  T nullImage( typename T::Domain( firstPoint, lastPoint ));
59 
60  HeaderField header[ MAX_HEADERNUMLINES ];
61 
62 #ifdef WIN32
63  errno_t err;
64  err = fopen_s( &fin, filename.c_str() , "rb" );
65  if ( err )
66  {
67  trace.error() << "VolReader : can't open " << filename << std::endl;
68  throw dgtalexception;
69  }
70 #else
71  fin = fopen( filename.c_str() , "rb" );
72 #endif
73 
74  if ( fin == NULL )
75  {
76  trace.error() << "VolReader : can't open " << filename << std::endl;
77  throw dgtalexception;
78  }
79 
80 
81  // Read header
82  // Buf for a line
83  char buf[128];
84  int linecount = 1;
85  int fieldcount = 0;
86 
87  // Read the file line by line until ".\n" is found
88  for ( char *line = fgets( buf, 128, fin );
89  line && strcmp( line, ".\n" ) != 0 ;
90  line = fgets( line, 128, fin ), ++linecount
91  )
92  {
93 
94  if ( line[strlen( line ) - 1] != '\n' )
95  {
96  trace.error() << "VolReader: Line " << linecount << " too long" << std::endl;
97  throw dgtalexception;
98  }
99 
100  int i;
101  for ( i = 0; line[i] && line[i] != ':'; ++i )
102  ;
103 
104  if ( i == 0 || i >= 126 || line[i] != ':' )
105  {
106  trace.error() << "VolReader: Invalid header read at line " << linecount << std::endl;
107  throw dgtalexception;
108  }
109  else
110  {
111 
112  if ( fieldcount == MAX_HEADERNUMLINES )
113  {
114  trace.warning() << "VolReader: Too many lines in HEADER, ignoring\n";
115  continue;
116  }
117  if ( fieldcount > MAX_HEADERNUMLINES )
118  continue;
119 
120  // Remove \n from end of line
121  if ( line[ strlen( line ) - 1 ] == '\n' )
122  line[ strlen( line ) - 1 ] = 0;
123 
124  // hack : split line into two str ...
125  line[i] = 0;
126  header[ fieldcount++ ] = HeaderField( line, line + i + 2 );
127  // +2 cause we skip the space
128  // following the colon
129  }
130  }
131 
132  // Check required headers
133  for ( int i = 0; requiredHeaders[i]; ++i )
134  {
135  if ( getHeaderValue( "Version" , header ) != NULL &&
136  ( strcmp( requiredHeaders[i], "Int-Endian" ) == 0 ||
137  strcmp( requiredHeaders[i], "Voxel-Endian" ) == 0 ) )
138  {
139  continue;
140  }
141  if ( getHeaderField( requiredHeaders[i] , header ) == -1 )
142  {
143  trace.error() << "VolReader: Required Header Field missing: "
144  << requiredHeaders[i] << std::endl;
145  throw dgtalexception;
146 
147  }
148  }
149 
150  int sx = 0, sy= 0, sz= 0;
151  int cx = 0, cy= 0, cz= 0;
152  int version = -1;
153 
154  getHeaderValueAsInt( "X", &sx, header );
155  getHeaderValueAsInt( "Y", &sy, header );
156  getHeaderValueAsInt( "Z", &sz, header );
157  getHeaderValueAsInt( "Version", &version, header);
158 
159  if (! ((version == 2) || (version == 3)))
160  {
161  trace.error() << "VolReader: invalid Version header (must be either 2 or 3)\n";
162  throw dgtalexception;
163  }
164 
165 
166  //Raw Data
167  if( getHeaderValueAsInt( "Center-X", &cx, header ) == 0 )
168  {
169  getHeaderValueAsInt( "Center-Y", &cy, header );
170  getHeaderValueAsInt( "Center-Z", &cz, header );
171 
172  firstPoint[0] = cx - (sx - 1)/2;
173  firstPoint[1] = cy - (sy - 1)/2;
174  firstPoint[2] = cz - (sz - 1)/2;
175  lastPoint[0] = cx + sx/2;
176  lastPoint[1] = cy + sy/2;
177  lastPoint[2] = cz + sz/2;
178  }
179  else
180  {
181  firstPoint = T::Point::zero;
182  lastPoint[0] = sx - 1;
183  lastPoint[1] = sy - 1;
184  lastPoint[2] = sz - 1;
185  }
186 
187  typename T::Domain domain( firstPoint, lastPoint );
188 
189  try
190  {
191  T image( domain );
192 
193  long count = 0;
194  unsigned char val;
195  typename T::Domain::ConstIterator it = domain.begin();
196  long int total = sx * sy * sz;
197  std::stringstream main;
198 
199  //main read loop
200  while (( count < total ) && ( fin ) )
201  {
202  val = getc( fin );
203  main << val;
204  count++;
205  }
206 
207  if ( count != total )
208  {
209  trace.error() << "VolReader: can't read file (raw data) !\n";
210  throw dgtalexception;
211  }
212 
213  //Uncompress if needed
214  if(version == 3)
215  {
216  std::stringstream uncompressed;
217  boost::iostreams::filtering_streambuf<boost::iostreams::input> in;
218  in.push(boost::iostreams::zlib_decompressor());
219  in.push(main);
220  boost::iostreams::copy(in, uncompressed);
221  //Apply to the image structure
222  for(auto i=0; i < total; ++i)
223  {
224  val = uncompressed.get();
225  image.setValue(( *it ), aFunctor(val) );
226  it++;
227  }
228  }
229  else
230  {
231  //Apply to the image structure
232  for(auto i=0; i < total; ++i)
233  {
234  val = main.get();
235  image.setValue(( *it ), aFunctor(val) );
236  it++;
237  }
238  }
239  fclose( fin );
240  return image;
241  }
242  catch ( ... )
243  {
244  trace.error() << "VolReader: not enough memory\n" ;
245  throw dgtalexception;
246  }
247 
248  }
249 
250 
251 
252  template <typename T, typename TFunctor>
253  const char *DGtal::VolReader<T, TFunctor>::requiredHeaders[] =
254  {
255  "X", "Y", "Z", "Voxel-Size", "Int-Endian", "Voxel-Endian", "Alpha-Color", NULL
256  };
257 
258 
259 
260 
261  template<typename T, typename TFunctor>
262  inline
263  int
264  DGtal::VolReader<T, TFunctor>::getHeaderField( const char *type, const HeaderField * header )
265  {
266 
267 
268  for ( int i = 0; i < MAX_HEADERNUMLINES; ++i )
269  {
270  if ( header[i].type != NULL && strcmp( header[i].type, type ) == 0 )
271  {
272  return i;
273  }
274  }
275  return -1;
276  }
277 
278 
279  template<typename T, typename TFunctor>
280  inline
281  const char *
282  DGtal::VolReader<T, TFunctor>::getHeaderValue( const char *type, const HeaderField * header )
283  {
284 
285  int i = getHeaderField( type, header );
286  if ( i == -1 )
287  return NULL;
288  return header[i].value;
289 
290  }
291 
292 
293  template<typename T, typename TFunctor>
294  inline
295  int
296  DGtal::VolReader<T, TFunctor>::getHeaderValueAsInt( const char *type, int *dest, const HeaderField * header )
297  {
298  int i = getHeaderField( type, header );
299  if ( i == -1 )
300  return 1;
301 
302  return sscanf( header[i].value, "%d", dest ) == 0;
303  }
304