4 #include <radarlib/radar.hpp>
8 using namespace OdimH5v21;
12 unsigned int_to_unsigned(
int val,
const char* desc)
17 errmsg +=
" has a negative value";
18 throw std::runtime_error(errmsg);
29 void ODIMLoader::request_quantity(
const std::string& name,
Scans<double>* volume)
31 to_load.insert(make_pair(name, volume));
34 void ODIMLoader::load(
const std::string& pathname)
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();
50 double range_scale = 0;
53 unsigned scan_count = int_to_unsigned(volume->getScanCount(),
"scan count");
54 double old_elevation = -1000.;
55 for (
unsigned src_elev = 0; src_elev < scan_count; ++src_elev)
57 unique_ptr<odim::PolarScan> scan(volume->getScan(src_elev));
58 double elevation = scan->getEAngle();
59 Available_Elevations.push_back(elevation);
63 if( elevation == old_elevation )
continue;
64 old_elevation=elevation;
67 range_scale = scan->getRangeScale();
69 double rs = scan->getRangeScale();
70 if (rs != range_scale)
72 LOG_ERROR(
"scan %d (elevation %f) has rangeScale %f that is different from %f in the previous scans",
73 src_elev, elevation, rs, range_scale);
74 throw runtime_error(
"rangeScale mismatch");
79 std::vector<odim::AZAngles> azangles = scan->getAzimuthAngles();
80 int rpm_sign = scan->getDirection();
82 std::vector<double> elevation_angles = scan->getElevationAngles();
84 unsigned beam_count = int_to_unsigned(scan->getNumRays(),
"number of rays");
85 if (azangles.size() != beam_count)
87 LOG_ERROR(
"elevation %f has %zd azimuth angles and %d rays", elevation, azangles.size(), beam_count);
88 throw runtime_error(
"mismatch between number of azumuth angles and number of rays");
91 unsigned beam_size = int_to_unsigned(scan->getNumBins(),
"beam size");
94 for (
auto& todo : to_load)
97 const string& name = todo.first;
101 if (!scan->hasQuantityData(name))
103 LOG_WARN(
"no %s found for elevation angle %f: skipping", name.c_str(), elevation);
108 unique_ptr<odim::PolarScanData> data(scan->getQuantityData(name));
112 target.
radarSite.height_r = volume->getAltitude()/1000.;
113 target.
radarSite.antennaTowerHeight = 0.;
115 vol_pol_scan.
nodata = data->getNodata() * data->getGain() + data->getOffset();
116 vol_pol_scan.
undetect = data->getUndetect() * data->getGain() + data->getOffset();
117 vol_pol_scan.
gain = data->getGain();
118 vol_pol_scan.
offset = data->getOffset();
119 vol_pol_scan.
cell_size = scan->getRangeScale();
122 odim::RayMatrix<double> matrix;
123 matrix.resize(beam_count, beam_size);
124 data->readTranslatedData(matrix);
126 for (
unsigned src_az = 0; src_az < beam_count; ++src_az)
128 Eigen::VectorXd beam(beam_size);
130 for (
unsigned i = 0; i < beam_size; ++i)
131 beam(i) = matrix.elem(src_az, i);
133 vol_pol_scan.row(src_az) = beam;
135 vol_pol_scan.
azimuths_real(src_az) = azangles[src_az].averagedAngle(rpm_sign);
139 for (
auto& i: to_load)
141 i.second->load_info = load_info;
145 void ODIMStorer::store(
const std::string& pathname)
147 namespace odim = OdimH5v21;
148 using namespace Radar;
151 shared_ptr<LoadInfo> load_info = make_shared<LoadInfo>();
152 load_info->filename = pathname;
154 unique_ptr<odim::OdimFactory> factory(
new odim::OdimFactory());
155 unique_ptr<odim::PolarVolume> volume(factory->openPolarVolume(pathname));
157 for(
unsigned i=0;i<to_store_int.size();i++)
158 for(
unsigned j=0;j<to_store_int[i]->size();j++)
160 vector<odim::PolarScan*> scans;
161 scans = volume->getScans(to_store_int[i]->scan(j).elevation,0.1);
162 shared_ptr<odim::PolarScan> scan;
163 if(scans.size()) scan.reset(scans[0]);
166 scan.reset(volume->createScan());
167 scan->setEAngle(to_store_int[i]->scan(j).elevation);
169 unique_ptr<odim::PolarScanData> data(scan->createQuantityData(to_store_int[i]->quantity));
170 data->setGain((
double)to_store_int[i]->scan(j).gain);
171 data->setOffset((
double)to_store_int[i]->scan(j).offset);
172 data->setNodata((
double)to_store_int[i]->scan(j).nodata);
173 data->setUndetect((
double)to_store_int[i]->scan(j).undetect);
174 odim::RayMatrix<unsigned short> matrix;
175 matrix.resize(to_store_int[i]->scan(j).beam_count,to_store_int[i]->scan(j).beam_size);
176 for(
unsigned ii=0;ii<to_store_int[i]->scan(j).beam_count;ii++)
177 for(
unsigned jj=0;jj<to_store_int[i]->scan(j).beam_size;jj++)
178 matrix.elem(ii,jj) = to_store_int[i]->scan(j)(ii,jj);
179 data->writeData(matrix);
181 for(
unsigned i=0;i<to_store_fp.size();i++)
182 for(
unsigned j=0;j<to_store_fp[i]->size();j++)
184 vector<odim::PolarScan*> scans;
185 scans = volume->getScans(to_store_fp[i]->scan(j).elevation,0.1);
186 shared_ptr<odim::PolarScan> scan;
187 if(scans.size()) scan.reset(scans[0]);
190 scan.reset(volume->createScan());
191 scan->setEAngle(to_store_fp[i]->scan(j).elevation);
193 unique_ptr<odim::PolarScanData> data(scan->createQuantityData(to_store_fp[i]->quantity));
194 data->setGain((
double)to_store_fp[i]->scan(j).gain);
195 data->setOffset((
double)to_store_fp[i]->scan(j).offset);
196 data->setNodata((
double)to_store_fp[i]->scan(j).nodata);
197 data->setUndetect((
double)to_store_fp[i]->scan(j).undetect);
198 odim::RayMatrix<float> matrix;
199 matrix.resize(to_store_fp[i]->scan(j).beam_count,to_store_fp[i]->scan(j).beam_size);
200 for(
unsigned ii=0;ii<to_store_fp[i]->scan(j).beam_count;ii++)
201 for(
unsigned jj=0;jj<to_store_fp[i]->scan(j).beam_size;jj++)
202 matrix.elem(ii,jj) = to_store_fp[i]->scan(j)(ii,jj);
203 data->writeAndTranslate(matrix,(
float)data->getOffset(),(float)data->getGain(),H5::PredType::NATIVE_UINT16);
206 for(
unsigned i=0;i<to_store_uchar.size();i++)
207 for(
unsigned j=0;j<to_store_uchar[i]->size();j++)
209 vector<odim::PolarScan*> scans;
210 scans = volume->getScans(to_store_uchar[i]->scan(j).elevation,0.1);
211 shared_ptr<odim::PolarScan> scan;
212 if(scans.size()) scan.reset(scans[0]);
215 scan.reset(volume->createScan());
216 scan->setEAngle(to_store_uchar[i]->scan(j).elevation);
218 unique_ptr<odim::PolarScanData> data(scan->createQuantityData(to_store_uchar[i]->quantity));
219 data->setGain((
double)to_store_uchar[i]->scan(j).gain);
220 data->setOffset((
double)to_store_uchar[i]->scan(j).offset);
221 data->setNodata((
double)to_store_uchar[i]->scan(j).nodata);
222 data->setUndetect((
double)to_store_uchar[i]->scan(j).undetect);
223 odim::RayMatrix<float> matrix;
224 matrix.resize(to_store_uchar[i]->scan(j).beam_count,to_store_uchar[i]->scan(j).beam_size);
225 for(
unsigned ii=0;ii<to_store_uchar[i]->scan(j).beam_count;ii++)
226 for(
unsigned jj=0;jj<to_store_uchar[i]->scan(j).beam_size;jj++)
227 matrix.elem(ii,jj) = to_store_uchar[i]->scan(j)(ii,jj);
228 data->writeAndTranslate(matrix,(
float)data->getOffset(),(float)data->getGain(),H5::PredType::NATIVE_UINT8);
231 for(
unsigned i=0;i<to_replace.size();i++)
232 for(
unsigned j=0;j<to_replace[i]->size();j++)
234 vector<odim::PolarScan*> scans;
235 scans = volume->getScans(to_replace[i]->scan(j).elevation,0.1);
237 shared_ptr<odim::PolarScan> scan;
238 if(scans.size()) scan.reset(scans[0]);
246 unique_ptr<odim::PolarScanData> data(scan->getQuantityData(to_replace[i]->quantity));
247 H5::AtomType type = data->getDataType();
249 odim::RayMatrix<double> matrix;
250 matrix.resize(to_replace[i]->scan(j).beam_count,to_replace[i]->scan(j).beam_size);
251 for(
unsigned ii=0;ii<to_replace[i]->scan(j).beam_count;ii++)
252 for(
unsigned jj=0;jj<to_replace[i]->scan(j).beam_size;jj++)
253 matrix.elem(ii,jj) = to_replace[i]->scan(j)(ii,jj);
254 data->writeAndTranslate(matrix,(
float)data->getOffset(),(float)data->getGain(),type);
260 void ODIMStorer::storeQuality(
const std::string& pathname,
const std::string & task,
bool RemoveQualityFields)
262 namespace odim = OdimH5v21;
263 using namespace Radar;
266 shared_ptr<LoadInfo> load_info = make_shared<LoadInfo>();
267 load_info->filename = pathname;
269 unique_ptr<odim::OdimFactory> factory(
new odim::OdimFactory());
270 unique_ptr<odim::PolarVolume> volume(factory->openPolarVolume(pathname));
272 for(
unsigned i=0;i<to_store_uchar.size();i++)
273 for(
unsigned j=0;j<to_store_uchar[i]->size();j++)
275 vector<odim::PolarScan*> scans;
276 scans = volume->getScans(to_store_uchar[i]->scan(j).elevation,0.1);
277 shared_ptr<odim::PolarScan> scan;
278 if(scans.size()) scan.reset(scans[0]);
281 scan.reset(volume->createScan());
282 scan->setEAngle(to_store_uchar[i]->scan(j).elevation);
284 if (RemoveQualityFields) {
285 int iqc =scan->getQualityCount();
286 for (
int iq=iqc-1; iq >=0 ; iq--){
287 scan->removeQuality(iq);
290 unique_ptr<odim::OdimQuality> quality(scan->createQuality());
291 quality->getWhat()->set(ATTRIBUTE_WHAT_OFFSET, 0.);
292 quality->getWhat()->set(ATTRIBUTE_WHAT_GAIN, 1.);
293 quality->getHow() ->set(ATTRIBUTE_HOW_TASK, task);
294 odim::RayMatrix<unsigned short> matrix;
295 matrix.resize(to_store_uchar[i]->scan(j).beam_count,to_store_uchar[i]->scan(j).beam_size);
296 for(
unsigned ii=0;ii<to_store_uchar[i]->scan(j).beam_count;ii++)
297 for(
unsigned jj=0;jj<to_store_uchar[i]->scan(j).beam_size;jj++)
298 matrix.elem(ii,jj) = to_store_uchar[i]->scan(j)(ii,jj);
299 quality->writeQuality(matrix);
double gain
Conversion factor.
double undetect
Minimum amount that can be measured.
Eigen::VectorXd elevations_real
Vector of actual elevations for each beam.
Codice per il caricamento di volumi ODIM in radarelab.
PolarScan - structure to describe a polarScan containing a matrix of data and conversion factors...
double cell_size
Size of a beam cell in meters.
double offset
Conversion factor.
PolarScan< T > & append_scan(unsigned beam_count, unsigned beam_size, double elevation, double cell_size)
Append a scan to this volume.
RadarSite radarSite
RadarSite.
std::string quantity
Odim quantity name.
Eigen::VectorXd azimuths_real
Vector of actual azimuths for each beam.
double nodata
Value used as 'no data' value.