4#include <radarlib/radar.hpp>
8using namespace OdimH5v21;
12unsigned int_to_unsigned(
int val,
const char* desc)
17 errmsg +=
" has a negative value";
18 throw std::runtime_error(errmsg);
31 to_load.insert(make_pair(name, volume));
36 LOG_CATEGORY(
"radar.io");
38 namespace odim = OdimH5v21;
39 using namespace Radar;
42 shared_ptr<LoadInfo> load_info = make_shared<LoadInfo>();
43 load_info->filename = pathname;
45 unique_ptr<odim::OdimFactory> factory(
new odim::OdimFactory());
46 unique_ptr<odim::PolarVolume> volume(factory->openPolarVolume(pathname));
48 load_info->acq_date = volume->getDateTime();
57 load_info->source_name = volume->getSource().Place.c_str() ;
59 catch (std::exception& stde)
61 std::cerr <<
"Errore durante l'esecuzione: " << stde.what() << std::endl;
65 std::cerr <<
"Errore sconosciuto" << std::endl;
68 double range_scale = 0;
71 unsigned scan_count = int_to_unsigned(volume->getScanCount(),
"scan count");
72 double old_elevation = -1000.;
73 for (
unsigned src_elev = 0; src_elev < scan_count; ++src_elev)
75 unique_ptr<odim::PolarScan> scan(volume->getScan(src_elev));
76 double elevation = scan->getEAngle();
77 Available_Elevations.push_back(elevation);
81 if( elevation == old_elevation )
continue;
82 old_elevation=elevation;
85 range_scale = scan->getRangeScale();
87 double rs = scan->getRangeScale();
88 if (rs != range_scale)
90 LOG_ERROR(
"scan %d (elevation %f) has rangeScale %f that is different from %f in the previous scans",
91 src_elev, elevation, rs, range_scale);
92 throw runtime_error(
"rangeScale mismatch");
97 std::vector<odim::AZAngles> azangles = scan->getAzimuthAngles();
98 int rpm_sign = scan->getDirection();
100 std::vector<double> elevation_angles = scan->getElevationAngles();
102 unsigned beam_count = int_to_unsigned(scan->getNumRays(),
"number of rays");
103 if (azangles.size() != beam_count)
105 LOG_ERROR(
"elevation %f has %zd azimuth angles and %d rays", elevation, azangles.size(), beam_count);
106 throw runtime_error(
"mismatch between number of azumuth angles and number of rays");
109 unsigned beam_size = int_to_unsigned(scan->getNumBins(),
"beam size");
115 const string& name = todo.first;
119 if (!scan->hasQuantityData(name))
121 LOG_WARN(
"no %s found for elevation angle %f: skipping", name.c_str(), elevation);
126 unique_ptr<odim::PolarScanData> data(scan->getQuantityData(name));
130 target.
radarSite.height_r = volume->getAltitude()/1000.;
131 target.
radarSite.antennaTowerHeight = 0.;
133 vol_pol_scan.
nodata = data->getNodata() * data->getGain() + data->getOffset();
134 vol_pol_scan.
undetect = data->getUndetect() * data->getGain() + data->getOffset();
135 vol_pol_scan.
gain = data->getGain();
136 vol_pol_scan.
offset = data->getOffset();
137 vol_pol_scan.
cell_size = scan->getRangeScale();
140 odim::RayMatrix<double> matrix;
141 matrix.resize(beam_count, beam_size);
142 data->readTranslatedData(matrix);
144 for (
unsigned src_az = 0; src_az < beam_count; ++src_az)
146 Eigen::VectorXd beam(beam_size);
148 for (
unsigned i = 0; i < beam_size; ++i)
149 beam(i) = matrix.elem(src_az, i);
151 vol_pol_scan.row(src_az) = beam;
153 vol_pol_scan.
azimuths_real(src_az) = azangles[src_az].averagedAngle(rpm_sign);
159 i.second->load_info = load_info;
163void ODIMStorer::store(
const std::string& pathname)
165 namespace odim = OdimH5v21;
166 using namespace Radar;
169 shared_ptr<LoadInfo> load_info = make_shared<LoadInfo>();
170 load_info->filename = pathname;
172 unique_ptr<odim::OdimFactory> factory(
new odim::OdimFactory());
173 unique_ptr<odim::PolarVolume> volume(factory->openPolarVolume(pathname));
175 for(
unsigned i=0;i<to_store_int.size();i++)
176 for(
unsigned j=0;j<to_store_int[i]->size();j++)
178 vector<odim::PolarScan*> scans;
179 scans = volume->getScans(to_store_int[i]->scan(j).elevation,0.1);
180 shared_ptr<odim::PolarScan> scan;
181 if(scans.size()) scan.reset(scans[0]);
184 scan.reset(volume->createScan());
185 scan->setEAngle(to_store_int[i]->scan(j).elevation);
187 unique_ptr<odim::PolarScanData> data(scan->createQuantityData(to_store_int[i]->quantity));
188 data->setGain((
double)to_store_int[i]->scan(j).gain);
189 data->setOffset((
double)to_store_int[i]->scan(j).offset);
190 data->setNodata((
double)to_store_int[i]->scan(j).nodata);
191 data->setUndetect((
double)to_store_int[i]->scan(j).undetect);
192 odim::RayMatrix<unsigned short> matrix;
193 matrix.resize(to_store_int[i]->scan(j).beam_count,to_store_int[i]->scan(j).beam_size);
194 for(
unsigned ii=0;ii<to_store_int[i]->scan(j).beam_count;ii++)
195 for(
unsigned jj=0;jj<to_store_int[i]->scan(j).beam_size;jj++)
196 matrix.elem(ii,jj) = to_store_int[i]->scan(j)(ii,jj);
197 data->writeData(matrix);
199 for(
unsigned i=0;i<to_store_fp.size();i++)
200 for(
unsigned j=0;j<to_store_fp[i]->size();j++)
202 vector<odim::PolarScan*> scans;
203 scans = volume->getScans(to_store_fp[i]->scan(j).elevation,0.1);
204 shared_ptr<odim::PolarScan> scan;
205 if(scans.size()) scan.reset(scans[0]);
208 scan.reset(volume->createScan());
209 scan->setEAngle(to_store_fp[i]->scan(j).elevation);
211 unique_ptr<odim::PolarScanData> data(scan->createQuantityData(to_store_fp[i]->quantity));
212 data->setGain((
double)to_store_fp[i]->scan(j).gain);
213 data->setOffset((
double)to_store_fp[i]->scan(j).offset);
214 data->setNodata((
double)to_store_fp[i]->scan(j).nodata);
215 data->setUndetect((
double)to_store_fp[i]->scan(j).undetect);
216 odim::RayMatrix<float> matrix;
217 matrix.resize(to_store_fp[i]->scan(j).beam_count,to_store_fp[i]->scan(j).beam_size);
218 for(
unsigned ii=0;ii<to_store_fp[i]->scan(j).beam_count;ii++)
219 for(
unsigned jj=0;jj<to_store_fp[i]->scan(j).beam_size;jj++)
220 matrix.elem(ii,jj) = to_store_fp[i]->scan(j)(ii,jj);
221 data->writeAndTranslate(matrix,(
float)data->getOffset(),(float)data->getGain(),H5::PredType::NATIVE_UINT16);
224 for(
unsigned i=0;i<to_store_uchar.size();i++)
225 for(
unsigned j=0;j<to_store_uchar[i]->size();j++)
227 vector<odim::PolarScan*> scans;
228 scans = volume->getScans(to_store_uchar[i]->scan(j).elevation,0.1);
229 shared_ptr<odim::PolarScan> scan;
230 if(scans.size()) scan.reset(scans[0]);
233 scan.reset(volume->createScan());
234 scan->setEAngle(to_store_uchar[i]->scan(j).elevation);
236 unique_ptr<odim::PolarScanData> data(scan->createQuantityData(to_store_uchar[i]->quantity));
237 data->setGain((
double)to_store_uchar[i]->scan(j).gain);
238 data->setOffset((
double)to_store_uchar[i]->scan(j).offset);
239 data->setNodata((
double)to_store_uchar[i]->scan(j).nodata);
240 data->setUndetect((
double)to_store_uchar[i]->scan(j).undetect);
241 odim::RayMatrix<float> matrix;
242 matrix.resize(to_store_uchar[i]->scan(j).beam_count,to_store_uchar[i]->scan(j).beam_size);
243 for(
unsigned ii=0;ii<to_store_uchar[i]->scan(j).beam_count;ii++)
244 for(
unsigned jj=0;jj<to_store_uchar[i]->scan(j).beam_size;jj++)
245 matrix.elem(ii,jj) = to_store_uchar[i]->scan(j)(ii,jj);
246 data->writeAndTranslate(matrix,(
float)data->getOffset(),(float)data->getGain(),H5::PredType::NATIVE_UINT8);
249 for(
unsigned i=0;i<to_replace.size();i++)
250 for(
unsigned j=0;j<to_replace[i]->size();j++)
252 vector<odim::PolarScan*> scans;
253 scans = volume->getScans(to_replace[i]->scan(j).elevation,0.1);
255 shared_ptr<odim::PolarScan> scan;
256 if(scans.size()) scan.reset(scans[0]);
264 unique_ptr<odim::PolarScanData> data(scan->getQuantityData(to_replace[i]->quantity));
265 H5::AtomType type = data->getDataType();
267 odim::RayMatrix<double> matrix;
268 matrix.resize(to_replace[i]->scan(j).beam_count,to_replace[i]->scan(j).beam_size);
269 for(
unsigned ii=0;ii<to_replace[i]->scan(j).beam_count;ii++)
270 for(
unsigned jj=0;jj<to_replace[i]->scan(j).beam_size;jj++)
271 matrix.elem(ii,jj) = to_replace[i]->scan(j)(ii,jj);
272 data->writeAndTranslate(matrix,(
float)data->getOffset(),(float)data->getGain(),type);
278void ODIMStorer::storeQuality(
const std::string& pathname,
const std::string & task,
bool RemoveQualityFields)
280 namespace odim = OdimH5v21;
281 using namespace Radar;
284 shared_ptr<LoadInfo> load_info = make_shared<LoadInfo>();
285 load_info->filename = pathname;
287 unique_ptr<odim::OdimFactory> factory(
new odim::OdimFactory());
288 unique_ptr<odim::PolarVolume> volume(factory->openPolarVolume(pathname));
290 for(
unsigned i=0;i<to_store_uchar.size();i++)
291 for(
unsigned j=0;j<to_store_uchar[i]->size();j++)
293 vector<odim::PolarScan*> scans;
294 scans = volume->getScans(to_store_uchar[i]->scan(j).elevation,0.1);
295 shared_ptr<odim::PolarScan> scan;
296 if(scans.size()) scan.reset(scans[0]);
299 scan.reset(volume->createScan());
300 scan->setEAngle(to_store_uchar[i]->scan(j).elevation);
302 if (RemoveQualityFields) {
303 int iqc =scan->getQualityCount();
304 for (
int iq=iqc-1; iq >=0 ; iq--){
305 scan->removeQuality(iq);
308 unique_ptr<odim::OdimQuality> quality(scan->createQuality());
309 quality->getWhat()->set(ATTRIBUTE_WHAT_OFFSET, 0.);
310 quality->getWhat()->set(ATTRIBUTE_WHAT_GAIN, 1.);
311 quality->getHow() ->set(ATTRIBUTE_HOW_TASK, task);
312 odim::RayMatrix<unsigned short> matrix;
313 matrix.resize(to_store_uchar[i]->scan(j).beam_count,to_store_uchar[i]->scan(j).beam_size);
314 for(
unsigned ii=0;ii<to_store_uchar[i]->scan(j).beam_count;ii++)
315 for(
unsigned jj=0;jj<to_store_uchar[i]->scan(j).beam_size;jj++)
316 matrix.elem(ii,jj) = to_store_uchar[i]->scan(j)(ii,jj);
317 quality->writeQuality(matrix);
double nodata
Value used as 'no data' value.
double undetect
Minimum amount that can be measured.
double offset
Conversion factor.
double gain
Conversion factor.
PolarScan - structure to describe a polarScan containing a matrix of data and conversion factors.
RadarSite radarSite
RadarSite.
std::string quantity
Odim quantity name.
PolarScan< T > & append_scan(unsigned beam_count, unsigned beam_size, double elevation, double cell_size, const T &default_value=algo::DBZ::BYTEtoDB(1))
Append a scan to this volume.
Sequence of PolarScans which can have a different beam count for each elevation.
Codice per il caricamento di volumi ODIM in radarelab.
Eigen::VectorXd elevations_real
Vector of actual elevations for each beam.
Eigen::VectorXd azimuths_real
Vector of actual azimuths for each beam.
double cell_size
Size of a beam cell in meters.
std::map< std::string, Scans< double > * > to_load
Map used to specify quantity to be loaded.
void request_quantity(const std::string &name, Scans< double > *volume)
Define a request - Fill to_load attribute
void load(const std::string &pathname)
Load method.