3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
rasterimagereader.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24 #ifndef DUMUX_RASTER_IMAGE_READER_HH
25 #define DUMUX_RASTER_IMAGE_READER_HH
26
27 #include <cassert>
28 #include <string>
29 #include <vector>
30 #include <fstream>
31 #include <sstream>
32 #include <algorithm>
33 #include <map>
34 #include <iterator>
35 #include <iostream>
36
37 #include <dune/common/exceptions.hh>
38
39namespace Dumux {
40
47{
48public:
52 struct Format
53 {
54 std::string magicNumber;
55 std::string type;
56 std::string encoding;
57 };
58
63 {
65 std::size_t nCols;
66 std::size_t nRows;
67 std::size_t maxValue = 1;
68 };
69
74 template<class T>
75 class Result : private std::vector<T>
76 {
77 using Parent = std::vector<T>;
78 public:
79 Result() = delete;
80
84 Result(const std::vector<T>& data, const HeaderData& header)
85 : Parent(data)
86 , header_(header)
87 {}
88
92 Result(std::vector<T>&& data, HeaderData&& header)
93 : Parent(std::move(data))
94 , header_(std::move(header))
95 {}
96
98 const HeaderData& header() const { return header_; }
99
100 // expose some methods of std::vector
101 using Parent::operator[];
102 using Parent::begin;
103 using Parent::end;
104 using Parent::size;
105
106 private:
107 HeaderData header_;
108 };
109
115 static Format getFormat(const std::string& magicNumber)
116 {
117 static const auto format = []{
118 std::map<std::string, Format> format;
119 format["P1"] = Format{"P1", "Portable BitMap", "ASCII"};
120 format["P2"] = Format{"P2", "Portable GrayMap", "ASCII"};
121 format["P3"] = Format{"P3", "Portable PixMap", "ASCII"};
122 format["P4"] = Format{"P4", "Portable BitMap", "binary"};
123 format["P5"] = Format{"P5", "Portable GrayMap", "binary"};
124 format["P5"] = Format{"P5", "Portable PixMap", "binary"};
125 return format;
126 }();
127
128 if (!format.count(magicNumber))
129 DUNE_THROW(Dune::IOError, magicNumber << " is not a valid magic number for the Netpbm format");
130
131 return format.at(magicNumber);
132 }
133
142 static Result<bool> readPBM(const std::string& fileName, const bool useDuneGridOrdering = true)
143 {
144 std::ifstream infile(fileName, std::ios::binary);
145
146 if (!infile)
147 DUNE_THROW(Dune::IOError, "Raster data file not found or corrupt");
148
149 HeaderData headerData = readHeader(infile);
150 std::vector<bool> values;
151
152 if (headerData.format.magicNumber == "P1")
153 values = readPBMDataASCII_(infile, headerData);
154 else if (headerData.format.magicNumber == "P4")
155 values = readPBMDataBinary_(infile, headerData);
156 else
157 DUNE_THROW(Dune::IOError, headerData.format.magicNumber << " not supported. Use readPBM for P1 and P4 or readPGM for P2 and P5");
158
159 Result<bool> result(std::move(values), std::move(headerData));
160 printInfo(result);
161
162 if (useDuneGridOrdering)
163 applyDuneGridOrdering(result);
164
165 return result;
166 }
167
180 template<class ValueType = std::uint8_t>
181 static Result<ValueType> readPGM(const std::string& fileName, const bool useDuneGridOrdering = true)
182 {
183 std::ifstream infile(fileName, std::ios::binary);
184
185 if (!infile)
186 DUNE_THROW(Dune::IOError, "Raster data file not found or corrupt");
187
188 HeaderData headerData = readHeader(infile);
189 std::vector<ValueType> values;
190
191 if (headerData.format.magicNumber == "P2")
192 values = NetPBMReader::template readPGMDataASCII_<ValueType>(infile, headerData);
193 else if (headerData.format.magicNumber == "P5")
194 values = NetPBMReader::template readPGMDataBinary_<ValueType>(infile, headerData);
195 else
196 DUNE_THROW(Dune::IOError, headerData.format.magicNumber << " not supported. Use readPBM for P1 and P4 or readPGM for P2 and P5");
197
198 Result<ValueType> result(std::move(values), std::move(headerData));
199 printInfo(result);
200
201 if (useDuneGridOrdering)
202 applyDuneGridOrdering(result);
203
204 return result;
205 }
206
212 static HeaderData readHeader(std::ifstream& infile)
213 {
214 HeaderData headerData;
215 std::string inputLine;
216
217 // First line : get format.
218 std::getline(infile, inputLine);
219 headerData.format = getFormat(inputLine);
220 const auto magicNumber = headerData.format.magicNumber;
221
222 // Read dimensions and maximum value (for non-b/w images).
223 while (!infile.eof())
224 {
225 std::getline(infile, inputLine);
226
227 auto isComment = [](const auto& s)
228 { return (s.find("#") != std::string::npos); };
229
230 // Skip comments.
231 if (isComment(inputLine))
232 continue;
233
234 // The first line after the comments contains the dimensions.
235 headerData.nCols = std::stoi(inputLine.substr(0, inputLine.find(" ")));
236 headerData.nRows = std::stoi(inputLine.substr(inputLine.find(" ") + 1));
237
238 // Grayscale images additionaly contain a maxium value in the header.
239 if (magicNumber != "P1" && magicNumber != "P4")
240 {
241 std::getline(infile, inputLine);
242 headerData.maxValue = std::stoi(inputLine);
243 }
244 break;
245 }
246 return headerData;
247 }
248
255 template<class T>
256 static void applyDuneGridOrdering(Result<T>& result)
257 {
258 auto tmp = result;
259 for (std::size_t i = 0; i < result.size(); i += result.header().nCols)
260 std::swap_ranges((result.begin() + i), (result.begin() + i + result.header().nCols), (tmp.end() - i - result.header().nCols));
261 }
262
268 template<class T>
269 static void printInfo(const Result<T>& result)
270 {
271 const Format& format = result.header().format;
272 std::cout << "Reading " << format.type << " file (" << format.encoding << ")" << std::endl;
273 std::cout << "Dimensions : " << result.header().nCols << " " << result.header().nRows << std::endl;
274 std::cout << "Maximum value : " << result.header().maxValue << std::endl;
275 }
276
284 template<class Image, class T>
285 static void fillImage(Image& image, const Result<T>& result)
286 {
287 const auto nCols = result.header().nCols;
288 const auto nRows = result.header().nRows;
289 using RowType = std::decay_t<decltype(image[0])>;
290 image.resize(nRows, RowType(nCols));
291
292 std::size_t rowIdx = 0;
293 std::size_t colIdx = 0;
294 for (const auto& val : result)
295 {
296 image[rowIdx][colIdx] = val;
297
298 // start a new line after nCols entries
299 if (++colIdx == nCols)
300 {
301 colIdx = 0;
302 ++rowIdx;
303 }
304 }
305 }
306
312 template<class Image>
313 static auto flattenImageToVector(const Image& image)
314 {
315 using OutputValueType = std::decay_t<decltype(image[0][0])>;
316
317 std::vector<OutputValueType> data;
318 data.reserve(image.size()*image[0].size());
319 for (const auto& row: image)
320 for (const auto& col : row)
321 data.push_back(col);
322
323 return data;
324 }
325
326private:
335 static std::vector<bool> readPBMDataASCII_(std::ifstream& infile,
336 const HeaderData& headerData)
337 {
338 std::string inputLine;
339 std::vector<bool> data;
340 data.reserve(numPixel_(headerData));
341
342 while (!infile.eof())
343 {
344 std::getline(infile, inputLine);
345 inputLine.erase(std::remove(inputLine.begin(), inputLine.end(), '\n'), inputLine.end());
346 inputLine.erase(std::remove(inputLine.begin(), inputLine.end(), ' '), inputLine.end());
347 if (!inputLine.empty())
348 {
349 for (const auto& value : inputLine)
350 {
351 assert(value == '0' || value == '1');
352 data.push_back(value - '0'); // convert char to int
353 }
354 }
355 }
356
357 return data;
358 }
359
368 static std::vector<bool> readPBMDataBinary_(std::ifstream& infile,
369 const HeaderData& headerData)
370 {
371 std::vector<bool> data(numPixel_(headerData));
372
373 std::size_t nBytes = 0;
374 std::size_t bitIndex = 0;
375 using Bit = std::uint8_t;
376 Bit b = 0;
377 for (std::size_t j = 0; j < headerData.nRows; j++)
378 {
379 for (std::size_t i = 0; i < headerData.nCols; i++)
380 {
381 if (i%8 == 0)
382 {
383 char tmp;
384 infile.read(&tmp, 1);
385 b = static_cast<Bit>(tmp);
386 if (infile.eof())
387 DUNE_THROW(Dune::IOError, "Failed reading byte " << nBytes);
388
389 ++nBytes;
390 }
391
392 const Bit k = 7 - (i % 8);
393 data[bitIndex++] = static_cast<bool>((b >> k) & 1);
394 }
395 }
396
397 return data;
398 }
399
412 template<class ValueType = std::uint8_t>
413 static std::vector<ValueType> readPGMDataASCII_(std::ifstream& infile,
414 const HeaderData& headerData)
415 {
416 std::string inputLine;
417
418 std::vector<ValueType> data;
419 data.reserve(numPixel_(headerData));
420
421 while (!infile.eof())
422 {
423 std::getline(infile, inputLine);
424 if (inputLine.empty())
425 continue;
426
427 // if the line contains multiple comma-separated values, store them individually in a vector
428 if (inputLine.find(" ") != std::string::npos)
429 {
430 std::istringstream iss(inputLine);
431 std::vector<std::string> tokens;
432 std::copy(std::istream_iterator<std::string>(iss),
433 std::istream_iterator<std::string>(),
434 std::back_inserter(tokens));
435
436 for (const auto& t : tokens)
437 data.push_back(std::stoi(t)); // convert string to integer type
438 }
439 else
440 data.push_back(std::stoi(inputLine));
441 }
442
443 return data;
444 }
445
458 template<class ValueType = std::uint8_t>
459 static std::vector<ValueType> readPGMDataBinary_(std::ifstream& infile,
460 const HeaderData& headerData)
461 {
462 std::string inputLine;
463
464 // check the size of the binary part of the file
465 const auto curPos = infile.tellg();
466 infile.seekg(0, std::ios::end);
467 const auto endPos = infile.tellg();
468 const auto size = endPos - curPos;
469 if (size != numPixel_(headerData))
470 DUNE_THROW(Dune::IOError, "Binary file size does not match with raster image size");
471
472 // reset to the current position
473 infile.seekg(curPos, std::ios::beg);
474
475 // extract the binary data
476 std::vector<std::uint8_t> data(size);
477 infile.read(reinterpret_cast<char*>(&data[0]), size);
478
479 // convert std::uint8_t to ValueType
480 return std::vector<ValueType>(data.begin(), data.end());
481 }
482
488 static std::size_t numPixel_(const HeaderData& headerData)
489 {
490 return headerData.nRows*headerData.nCols;
491 }
492};
493
494} // namespace Dumux
495
496#endif // DUMUX_RASTER_IMAGE_READER_HH
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
A simple reader class for the Netpbm format (https://en.wikipedia.org/wiki/Netpbm_format)....
Definition: rasterimagereader.hh:47
static HeaderData readHeader(std::ifstream &infile)
Returns the header data of the image file.
Definition: rasterimagereader.hh:212
static void fillImage(Image &image, const Result< T > &result)
Fill a pre-defined 2D image object, e.g. std::vector<std::vector<bool>>, with the pixel values stored...
Definition: rasterimagereader.hh:285
static Format getFormat(const std::string &magicNumber)
A helper function to retrieve the format data from a given magic number.
Definition: rasterimagereader.hh:115
static Result< bool > readPBM(const std::string &fileName, const bool useDuneGridOrdering=true)
Reads a *pbm (black and white) in ASCII or binary encoding. Returns a struct that contains both the p...
Definition: rasterimagereader.hh:142
static auto flattenImageToVector(const Image &image)
Flattens a 2D image object to a 1D container.
Definition: rasterimagereader.hh:313
static void applyDuneGridOrdering(Result< T > &result)
Change the ordering of the pixels according to Dune's convention, shifting the origin from upper left...
Definition: rasterimagereader.hh:256
static Result< ValueType > readPGM(const std::string &fileName, const bool useDuneGridOrdering=true)
Reads a *.pgm (grayscale) in ASCII or binary encoding. Returns a struct that contains both the pixel ...
Definition: rasterimagereader.hh:181
static void printInfo(const Result< T > &result)
Print the data contained in the header.
Definition: rasterimagereader.hh:269
A struct that holds all information of the image format.
Definition: rasterimagereader.hh:53
std::string encoding
Definition: rasterimagereader.hh:56
std::string magicNumber
Definition: rasterimagereader.hh:54
std::string type
Definition: rasterimagereader.hh:55
A struct that contains all header data of the image.
Definition: rasterimagereader.hh:63
Format format
Definition: rasterimagereader.hh:64
std::size_t maxValue
Definition: rasterimagereader.hh:67
std::size_t nCols
Definition: rasterimagereader.hh:65
std::size_t nRows
Definition: rasterimagereader.hh:66
The return type of the reading functions. Holds the actual pixel values and the header data.
Definition: rasterimagereader.hh:76
const HeaderData & header() const
Returns the header data.
Definition: rasterimagereader.hh:98
Result(const std::vector< T > &data, const HeaderData &header)
Contruct from data and header by copy.
Definition: rasterimagereader.hh:84
Result(std::vector< T > &&data, HeaderData &&header)
Contruct from data and header by move.
Definition: rasterimagereader.hh:92