OGS
anonymous_namespace{NetCDFRasterReader.cpp} Namespace Reference

Functions

GeoLib::RasterHeader readHeaderFromNetCDF (std::multimap< std::string, netCDF::NcVar > const &variables)
 
std::vector< double > readDataFromNetCDF (std::multimap< std::string, netCDF::NcVar > const &variables, std::string_view const var_name, std::size_t const dimension_number, GeoLib::MinMaxPoints const &min_max_points, GeoLib::RasterHeader &header)
 
std::unique_ptr< GeoLib::RasterreadNetCDF (std::filesystem::path const &filepath, std::string_view const var_name, std::size_t const dimension_number, GeoLib::MinMaxPoints const &min_max_points)
 
GeoLib::NamedRaster readRasterFromFile (std::filesystem::path const &path, std::filesystem::path filename, std::string const &var_name, std::size_t const dimension_number, GeoLib::MinMaxPoints const &min_max_points)
 

Function Documentation

◆ readDataFromNetCDF()

std::vector< double > anonymous_namespace{NetCDFRasterReader.cpp}::readDataFromNetCDF ( std::multimap< std::string, netCDF::NcVar > const & variables,
std::string_view const var_name,
std::size_t const dimension_number,
GeoLib::MinMaxPoints const & min_max_points,
GeoLib::RasterHeader & header )

Definition at line 74 of file NetCDFRasterReader.cpp.

78{
79 // second read raster header values
80 std::vector<double> raster_data;
81 for (auto [variable_name, variable] : variables)
82 {
83 // uncomment for debugging
84 // DBUG("checking variable '{}'", variable_name);
85 if (variable_name == var_name)
86 {
87 auto const n_dims = variable.getDimCount();
88 if (n_dims == 0)
89 {
90 continue;
91 }
92 // DBUG("variable '{}' has {} dimensions", variable_name, n_dims);
93
94 std::vector<std::size_t> counts(n_dims, 1);
95 for (int i = 0; i < n_dims; ++i)
96 {
97 counts[i] = variable.getDim(i).getSize();
98 }
99 // DBUG("counts {}", fmt::join(counts, ", "));
100 std::vector<std::size_t> start(n_dims, 0);
101 if (dimension_number > counts[0])
102 {
103 OGS_FATAL(
104 "variable '{}' has {} dimensions, requested dimension "
105 "number is {}",
106 variable_name, counts[0], dimension_number);
107 }
108 // time dimension
109 start[0] = dimension_number; // number of time step
110 counts[0] = 1; // read one time step
111
112 // y dimension - netcdf dimension number 1
113 // With the counts[1] and cell_size the 'usual' lower y-coordinate
114 // of the origin is computed.
115 // DBUG("reset header.origin[1]: original y0: {}, new y0: {} ",
116 // header.origin[1],
117 // header.origin[1] - counts[1] * header.cell_size);
118 header.origin[1] -= counts[1] * header.cell_size;
119
120 counts[2] = static_cast<int>(
121 std::ceil((min_max_points.max[0] - min_max_points.min[0]) /
122 header.cell_size) +
123 1);
124 counts[1] = static_cast<int>(
125 std::ceil((min_max_points.max[1] - min_max_points.min[1]) /
126 header.cell_size) +
127 1);
128 // x dimension - netcdf dimension number 2
129 start[2] = static_cast<int>(std::floor(
130 (min_max_points.min[0] - header.origin[0]) / header.cell_size));
131 // y dimension - netcdf dimension number 1
132 start[1] = static_cast<int>(std::floor(
133 (min_max_points.min[1] - header.origin[1]) / header.cell_size));
134
135 // DBUG("cut-out raster: pixel in x dimension: {} ", counts[2]);
136 // DBUG("cut-out raster: pixel in y dimension: {} ", counts[1]);
137 // DBUG("cut-out raster: start index x dimension: {}", start[2]);
138 // DBUG("cut-out raster: start index y dimension: {}", start[1]);
139
140 header.n_cols = counts[2];
141 header.n_rows = counts[1];
142 header.origin[0] += start[2] * header.cell_size;
143 // reset header y origin for cut out raster
144 header.origin[1] += start[1] * header.cell_size;
145 std::size_t const total_length =
146 std::accumulate(counts.cbegin(), counts.cend(), 1,
147 std::multiplies<std::size_t>());
148 raster_data.resize(total_length);
149 variable.getVar(start, counts, raster_data.data());
150
151 std::replace(raster_data.begin(), raster_data.end(), header.no_data,
152 0.0);
153 }
154 }
155 return raster_data;
156}
#define OGS_FATAL(...)
Definition Error.h:26
MathLib::Point3d origin
Definition Raster.h:32
std::size_t n_cols
Definition Raster.h:29
std::size_t n_rows
Definition Raster.h:30

References GeoLib::RasterHeader::cell_size, GeoLib::MinMaxPoints::max, GeoLib::MinMaxPoints::min, GeoLib::RasterHeader::n_cols, GeoLib::RasterHeader::n_rows, GeoLib::RasterHeader::no_data, OGS_FATAL, and GeoLib::RasterHeader::origin.

Referenced by readNetCDF().

◆ readHeaderFromNetCDF()

GeoLib::RasterHeader anonymous_namespace{NetCDFRasterReader.cpp}::readHeaderFromNetCDF ( std::multimap< std::string, netCDF::NcVar > const & variables)

Definition at line 27 of file NetCDFRasterReader.cpp.

29{
31
32 // remark: netcdf raster starts in the left upper corner
33 double cell_size_y = 0;
34 // first read:
35 // - origin, cell_size from GeoTransform attribute
36 // - header.no_data from _FillValue attribute
37 for (auto [variable_name, variable] : variables)
38 {
39 if (variable.isNull())
40 {
41 OGS_FATAL("Variable '{}' not found in file.", variable_name);
42 }
43 auto const& attributes = variable.getAtts();
44 for (auto [name, attribute] : attributes)
45 {
46 if (name == "GeoTransform")
47 {
48 std::vector<char> attribute_values;
49 attribute_values.resize(attribute.getAttLength());
50 attribute.getValues(attribute_values.data());
51 std::string const str_attribute_values(attribute_values.begin(),
52 attribute_values.end());
53 auto const strings = BaseLib::splitString(str_attribute_values);
54 header.origin[0] = std::stod(strings[0]); // origin x
55 header.cell_size = std::stod(strings[1]); // cell size x
56 header.origin[1] = std::stod(strings[3]); // origin y
57 cell_size_y = std::stod(strings[5]); // cell size y
58 DBUG(
59 "GeoTransform attribute values: origin_x: {}, origin_y: "
60 "{}, cell_size_x: {}, cell_size_y: {}",
61 header.origin[0], header.origin[1], header.cell_size,
62 cell_size_y);
63 }
64 if (name == "_FillValue")
65 {
66 attribute.getValues(&header.no_data);
67 }
68 }
69 }
70
71 return header;
72}
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
std::vector< std::string > splitString(std::string const &str)
Contains the relevant information when storing a geoscientific raster data.
Definition Raster.h:28

References GeoLib::RasterHeader::cell_size, DBUG(), GeoLib::RasterHeader::no_data, OGS_FATAL, GeoLib::RasterHeader::origin, and BaseLib::splitString().

Referenced by readNetCDF().

◆ readNetCDF()

std::unique_ptr< GeoLib::Raster > anonymous_namespace{NetCDFRasterReader.cpp}::readNetCDF ( std::filesystem::path const & filepath,
std::string_view const var_name,
std::size_t const dimension_number,
GeoLib::MinMaxPoints const & min_max_points )

Definition at line 158 of file NetCDFRasterReader.cpp.

163{
164 netCDF::NcFile dataset(filepath.string(), netCDF::NcFile::read);
165 if (dataset.isNull())
166 {
167 OGS_FATAL("Error opening file '{}'.", filepath.string());
168 }
169
170 auto const& variables = dataset.getVars();
171 GeoLib::RasterHeader header = readHeaderFromNetCDF(variables);
172 std::vector<double> raster_data = readDataFromNetCDF(
173 variables, var_name, dimension_number, min_max_points, header);
174
175 return std::make_unique<GeoLib::Raster>(header, raster_data.begin(),
176 raster_data.end());
177}
std::vector< double > readDataFromNetCDF(std::multimap< std::string, netCDF::NcVar > const &variables, std::string_view const var_name, std::size_t const dimension_number, GeoLib::MinMaxPoints const &min_max_points, GeoLib::RasterHeader &header)
GeoLib::RasterHeader readHeaderFromNetCDF(std::multimap< std::string, netCDF::NcVar > const &variables)

References OGS_FATAL, readDataFromNetCDF(), and readHeaderFromNetCDF().

Referenced by readRasterFromFile().

◆ readRasterFromFile()

GeoLib::NamedRaster anonymous_namespace{NetCDFRasterReader.cpp}::readRasterFromFile ( std::filesystem::path const & path,
std::filesystem::path filename,
std::string const & var_name,
std::size_t const dimension_number,
GeoLib::MinMaxPoints const & min_max_points )

Definition at line 180 of file NetCDFRasterReader.cpp.

186{
187 INFO("readRasterFromFile: '{}/{}'", path.string(), filename.string());
188
189 if (filename.extension() == ".nc")
190 {
191#ifdef OGS_USE_NETCDF
192 auto raster = readNetCDF(path / filename, var_name, dimension_number,
193 min_max_points);
194
195 return GeoLib::NamedRaster{filename.replace_extension().string() + "_" +
196 var_name + "_" +
197 std::to_string(dimension_number),
198 std::move(raster)};
199#else
200 OGS_FATAL("OGS was not build with NetCDF support. Can not read {}",
201 (path / filename).string());
202#endif
203 }
204 auto raster = std::unique_ptr<GeoLib::Raster>(
205 FileIO::AsciiRasterInterface::readRaster((path / filename).string()));
206 if (raster == nullptr)
207 {
208 OGS_FATAL("Could not read raster from file '{}'.",
209 (path / filename).string());
210 }
211 return GeoLib::NamedRaster{filename.replace_extension().string() + "_" +
212 var_name + "_" +
213 std::to_string(dimension_number),
214 std::move(raster)};
215}
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
static GeoLib::Raster * readRaster(std::string const &fname)
std::unique_ptr< GeoLib::Raster > readNetCDF(std::filesystem::path const &filepath, std::string_view const var_name, std::size_t const dimension_number, GeoLib::MinMaxPoints const &min_max_points)

References INFO(), OGS_FATAL, readNetCDF(), and FileIO::AsciiRasterInterface::readRaster().