3.5-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>
39
40namespace Dumux {
41
48{
49public:
53 struct Format
54 {
55 std::string magicNumber;
56 std::string type;
57 std::string encoding;
58 };
59
64 {
66 std::size_t nCols;
67 std::size_t nRows;
68 std::size_t maxValue = 1;
69 };
70
75 template<class T>
76 class Result : private std::vector<T>
77 {
78 using Parent = std::vector<T>;
79 public:
80 Result() = delete;
81
85 Result(const std::vector<T>& data, const HeaderData& header)
86 : Parent(data)
87 , header_(header)
88 {}
89
93 Result(std::vector<T>&& data, HeaderData&& header)
94 : Parent(std::move(data))
95 , header_(std::move(header))
96 {}
97
99 const HeaderData& header() const { return header_; }
100
101 // expose some methods of std::vector
102 using Parent::operator[];
103 using Parent::begin;
104 using Parent::end;
105 using Parent::size;
106
107 private:
108 HeaderData header_;
109 };
110
116 static Format getFormat(const std::vector<std::string_view>& firstLineTokes)
117 {
118 static const auto format = []{
119 std::map<std::string, Format> format;
120 format["P1"] = Format{"P1", "Portable BitMap", "ASCII"};
121 format["P2"] = Format{"P2", "Portable GrayMap", "ASCII"};
122 format["P3"] = Format{"P3", "Portable PixMap", "ASCII"};
123 format["P4"] = Format{"P4", "Portable BitMap", "binary"};
124 format["P5"] = Format{"P5", "Portable GrayMap", "binary"};
125 format["P6"] = Format{"P6", "Portable PixMap", "binary"};
126 return format;
127 }();
128
129 const std::string& magicNumber = std::string(firstLineTokes[0]);
130
131 if (!format.count(magicNumber))
132 DUNE_THROW(Dune::IOError, magicNumber << " is not a valid magic number for the Netpbm format");
133
134 return format.at(magicNumber);
135 }
136
145 static Result<bool> readPBM(const std::string& fileName, const bool useDuneGridOrdering = true)
146 {
147 std::ifstream infile(fileName, std::ios::binary);
148
149 if (!infile)
150 DUNE_THROW(Dune::IOError, "Raster data file not found or corrupt");
151
152 HeaderData headerData = readHeader(infile);
153 std::vector<bool> values;
154
155 if (headerData.format.magicNumber == "P1")
156 values = readPBMDataASCII_(infile, headerData);
157 else if (headerData.format.magicNumber == "P4")
158 values = readPBMDataBinary_(infile, headerData);
159 else
160 DUNE_THROW(Dune::IOError, headerData.format.magicNumber << " not supported. Use readPBM for P1 and P4 or readPGM for P2 and P5");
161
162 Result<bool> result(std::move(values), std::move(headerData));
163 printInfo(result);
164
165 if (useDuneGridOrdering)
166 applyDuneGridOrdering(result);
167
168 return result;
169 }
170
183 template<class ValueType = std::uint8_t>
184 static Result<ValueType> readPGM(const std::string& fileName, const bool useDuneGridOrdering = true)
185 {
186 std::ifstream infile(fileName, std::ios::binary);
187
188 if (!infile)
189 DUNE_THROW(Dune::IOError, "Raster data file not found or corrupt");
190
191 HeaderData headerData = readHeader(infile);
192 std::vector<ValueType> values;
193
194 if (headerData.format.magicNumber == "P2")
195 values = NetPBMReader::template readPGMDataASCII_<ValueType>(infile, headerData);
196 else if (headerData.format.magicNumber == "P5")
197 values = NetPBMReader::template readPGMDataBinary_<ValueType>(infile, headerData);
198 else
199 DUNE_THROW(Dune::IOError, headerData.format.magicNumber << " not supported. Use readPBM for P1 and P4 or readPGM for P2 and P5");
200
201 Result<ValueType> result(std::move(values), std::move(headerData));
202 printInfo(result);
203
204 if (useDuneGridOrdering)
205 applyDuneGridOrdering(result);
206
207 return result;
208 }
209
215 static HeaderData readHeader(std::ifstream& infile)
216 {
217 HeaderData headerData;
218 std::string inputLine;
219 std::size_t lineNumber = 0;
220
221 // First line : get format.
222 std::getline(infile, inputLine);
223 ++lineNumber;
224
225 const auto firstLineTokens = tokenize(inputLine, " ");
226 headerData.format = getFormat(firstLineTokens);
227 const auto magicNumber = headerData.format.magicNumber;
228
229 // dimensions could be given right after magic number (in same line)
230 if (firstLineTokens.size() > 2)
231 {
232 if (isBlackAndWhite_(magicNumber) && firstLineTokens.size() != 3)
233 DUNE_THROW(Dune::IOError, "Could not read first line for B/W image");
234
235 headerData.nCols = std::stoi(std::string(firstLineTokens[1]));
236 headerData.nRows = std::stoi(std::string(firstLineTokens[2]));
237
238 if (isGrayScale_(magicNumber))
239 {
240 if (firstLineTokens.size() == 4)
241 headerData.maxValue = std::stoi(std::string(firstLineTokens[3]));
242 if (firstLineTokens.size() > 4)
243 DUNE_THROW(Dune::IOError, "Could not read first line for grayscale image");
244 }
245 }
246 else
247 {
248 // Read dimensions and maximum value (for non-b/w images).
249 while (!infile.eof())
250 {
251 std::getline(infile, inputLine);
252 ++lineNumber;
253
254 // Skip comments.
255 if (isComment_(inputLine))
256 continue;
257
258 const auto tokens = tokenize(inputLine, " ");
259
260 // The first line after the comments contains the dimensions.
261 if (tokens.size() != 2)
262 DUNE_THROW(Dune::IOError, "Expecting " << [](auto size){ return size < 2 ? "both" : "only"; }(tokens.size()) << " dimensions (2 numbers) in line " << lineNumber);
263
264 headerData.nCols = std::stoi(std::string(tokens[0]));
265 headerData.nRows = std::stoi(std::string(tokens[1]));
266
267 // Grayscale images additionaly contain a maxium value in the header.
268 if (isGrayScale_(magicNumber))
269 {
270 std::getline(infile, inputLine);
271 ++lineNumber;
272
273 const auto token = tokenize(inputLine, " ");
274 if (token.size() != 1)
275 DUNE_THROW(Dune::IOError, "Expecting" << [](auto size){ return size == 0 ? "" : " only"; }(token.size()) << " intensity (one number) in line " << lineNumber);
276
277 headerData.maxValue = std::stoi(std::string(token[0]));
278 }
279 break;
280 }
281 }
282
283 return headerData;
284 }
285
292 template<class T>
293 static void applyDuneGridOrdering(Result<T>& result)
294 {
295 auto tmp = result;
296 for (std::size_t i = 0; i < result.size(); i += result.header().nCols)
297 std::swap_ranges((result.begin() + i), (result.begin() + i + result.header().nCols), (tmp.end() - i - result.header().nCols));
298 }
299
305 template<class T>
306 static void printInfo(const Result<T>& result)
307 {
308 const Format& format = result.header().format;
309 std::cout << "Reading " << format.type << " file (" << format.encoding << ")" << std::endl;
310 std::cout << "Dimensions : " << result.header().nCols << " " << result.header().nRows << std::endl;
311 std::cout << "Maximum value : " << result.header().maxValue << std::endl;
312 }
313
321 template<class Image, class T>
322 static void fillImage(Image& image, const Result<T>& result)
323 {
324 const auto nCols = result.header().nCols;
325 const auto nRows = result.header().nRows;
326 using RowType = std::decay_t<decltype(image[0])>;
327 image.resize(nRows, RowType(nCols));
328
329 std::size_t rowIdx = 0;
330 std::size_t colIdx = 0;
331 for (const auto val : result)
332 {
333 image[rowIdx][colIdx] = val;
334
335 // start a new line after nCols entries
336 if (++colIdx == nCols)
337 {
338 colIdx = 0;
339 ++rowIdx;
340 }
341 }
342 }
343
349 template<class Image>
350 static auto flattenImageToVector(const Image& image)
351 {
352 using OutputValueType = std::decay_t<decltype(image[0][0])>;
353
354 std::vector<OutputValueType> data;
355 data.reserve(image.size()*image[0].size());
356 for (const auto& row: image)
357 for (const auto col : row)
358 data.push_back(col);
359
360 return data;
361 }
362
363private:
364
365 static bool isBlackAndWhite_(const std::string& magicNumber)
366 {
367 return magicNumber == "P1" || magicNumber == "P4";
368 }
369
370 static bool isGrayScale_(const std::string& magicNumber)
371 {
372 return magicNumber == "P2" || magicNumber == "P5";
373 }
374
383 static std::vector<bool> readPBMDataASCII_(std::ifstream& infile,
384 const HeaderData& headerData)
385 {
386 std::string inputLine;
387 std::vector<bool> data;
388 data.reserve(numPixel_(headerData));
389
390 while (!infile.eof())
391 {
392 std::getline(infile, inputLine);
393 if (!isComment_(inputLine))
394 {
395 inputLine.erase(std::remove_if(inputLine.begin(), inputLine.end(), [](unsigned char c){ return std::isspace(c); }), inputLine.end());
396 if (!inputLine.empty())
397 {
398 for (const auto& value : inputLine)
399 {
400 assert(value == '0' || value == '1');
401 data.push_back(value - '0'); // convert char to int
402 }
403 }
404 }
405 }
406
407 return data;
408 }
409
418 static std::vector<bool> readPBMDataBinary_(std::ifstream& infile,
419 const HeaderData& headerData)
420 {
421 std::vector<bool> data(numPixel_(headerData));
422
423 // Skip potentially remaining comments in header section
424 // before reading binary content. We detect a comment by
425 // reading a line with std::getline and checking the resulting string.
426 // We continue reading new lines until no more comments are found. Then we
427 // need to set infile's current position to one line before the actual binary
428 // content, otherwise the following steps will fail.
429 std::string inputLine;
430 while (!infile.eof())
431 {
432 // store position before calling std::getline
433 const auto lastPos = infile.tellg();
434 std::getline(infile, inputLine);
435
436 // stop the loop if no more comment is found and go back one line
437 if (!isComment_(inputLine))
438 {
439 infile.seekg(lastPos);
440 break;
441 }
442 }
443
444 // read actual binary content
445 std::size_t nBytes = 0;
446 std::size_t bitIndex = 0;
447 using Bit = std::uint8_t;
448 Bit b = 0;
449 for (std::size_t j = 0; j < headerData.nRows; j++)
450 {
451 for (std::size_t i = 0; i < headerData.nCols; i++)
452 {
453 if (i%8 == 0)
454 {
455 char tmp;
456 infile.read(&tmp, 1);
457 b = static_cast<Bit>(tmp);
458 if (infile.eof())
459 DUNE_THROW(Dune::IOError, "Failed reading byte " << nBytes);
460
461 ++nBytes;
462 }
463
464 const Bit k = 7 - (i % 8);
465 data[bitIndex++] = static_cast<bool>((b >> k) & 1);
466 }
467 }
468
469 return data;
470 }
471
484 template<class ValueType = std::uint8_t>
485 static std::vector<ValueType> readPGMDataASCII_(std::ifstream& infile,
486 const HeaderData& headerData)
487 {
488 std::string inputLine;
489
490 std::vector<ValueType> data;
491 data.reserve(numPixel_(headerData));
492
493 while (!infile.eof())
494 {
495 std::getline(infile, inputLine);
496 if (inputLine.empty())
497 continue;
498
499 // if the line contains multiple comma-separated values, store them individually in a vector
500 if (inputLine.find(" ") != std::string::npos)
501 {
502 std::istringstream iss(inputLine);
503 std::vector<std::string> tokens;
504 std::copy(std::istream_iterator<std::string>(iss),
505 std::istream_iterator<std::string>(),
506 std::back_inserter(tokens));
507
508 for (const auto& t : tokens)
509 data.push_back(std::stoi(t)); // convert string to integer type
510 }
511 else
512 data.push_back(std::stoi(inputLine));
513 }
514
515 return data;
516 }
517
530 template<class ValueType = std::uint8_t>
531 static std::vector<ValueType> readPGMDataBinary_(std::ifstream& infile,
532 const HeaderData& headerData)
533 {
534 std::string inputLine;
535
536 // check the size of the binary part of the file
537 const auto curPos = infile.tellg();
538 infile.seekg(0, std::ios::end);
539 const auto endPos = infile.tellg();
540 const auto size = endPos - curPos;
541 if (size != numPixel_(headerData))
542 DUNE_THROW(Dune::IOError, "Binary file size does not match with raster image size");
543
544 // reset to the current position
545 infile.seekg(curPos, std::ios::beg);
546
547 // extract the binary data
548 std::vector<std::uint8_t> data(size);
549 infile.read(reinterpret_cast<char*>(&data[0]), size);
550
551 // convert std::uint8_t to ValueType
552 return std::vector<ValueType>(data.begin(), data.end());
553 }
554
560 static std::size_t numPixel_(const HeaderData& headerData)
561 {
562 return headerData.nRows*headerData.nCols;
563 }
564
568 static bool isComment_(const std::string_view line)
569 {
570 return line[0] == '#';
571 }
572};
573
574} // namespace Dumux
575
576#endif
Helpers for working with strings.
Definition: adapt.hh:29
std::vector< std::string_view > tokenize(std::string_view str, std::string_view delim)
Definition: stringutilities.hh:50
constexpr Line line
Definition: couplingmanager1d3d_line.hh:43
A simple reader class for the Netpbm format (https://en.wikipedia.org/wiki/Netpbm_format)....
Definition: rasterimagereader.hh:48
static HeaderData readHeader(std::ifstream &infile)
Returns the header data of the image file.
Definition: rasterimagereader.hh:215
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:322
static Format getFormat(const std::vector< std::string_view > &firstLineTokes)
A helper function to retrieve the format from tokens of the file's first line.
Definition: rasterimagereader.hh:116
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:145
static auto flattenImageToVector(const Image &image)
Flattens a 2D image object to a 1D container.
Definition: rasterimagereader.hh:350
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:293
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:184
static void printInfo(const Result< T > &result)
Print the data contained in the header.
Definition: rasterimagereader.hh:306
A struct that holds all information of the image format.
Definition: rasterimagereader.hh:54
std::string encoding
Definition: rasterimagereader.hh:57
std::string magicNumber
Definition: rasterimagereader.hh:55
std::string type
Definition: rasterimagereader.hh:56
A struct that contains all header data of the image.
Definition: rasterimagereader.hh:64
Format format
Definition: rasterimagereader.hh:65
std::size_t maxValue
Definition: rasterimagereader.hh:68
std::size_t nCols
Definition: rasterimagereader.hh:66
std::size_t nRows
Definition: rasterimagereader.hh:67
The return type of the reading functions. Holds the actual pixel values and the header data.
Definition: rasterimagereader.hh:77
const HeaderData & header() const
Returns the header data.
Definition: rasterimagereader.hh:99
Result(const std::vector< T > &data, const HeaderData &header)
Contruct from data and header by copy.
Definition: rasterimagereader.hh:85
Result(std::vector< T > &&data, HeaderData &&header)
Contruct from data and header by move.
Definition: rasterimagereader.hh:93