Elaboradar  0.1
 Tutto Classi Namespace File Funzioni Variabili Tipi enumerati (enum) Gruppi
elabora_volume.h
Vai alla documentazione di questo file.
1 
5 #include <radarelab/volume.h>
7 
8 using namespace std;
9 
10 namespace radarelab {
11 
12 namespace {bool check_undetect=false;}
13 
22 template<typename T>
23 inline bool good(const PolarScan<T>& scan, unsigned row, unsigned col)
24 {
25  if(check_undetect)
26  if(scan(row,col)!=scan.nodata&&scan(row,col)!=scan.undetect)
27  return true;
28  else return false;
29  else
30  if(scan(row,col)!=scan.nodata)
31  return true;
32  else return false;
33 }
34 
35 //===========================PolarScan manipulators==========================//
36 
43 template<typename T>
44 PolarScan<T> make_slope_scan(const PolarScan<T>& raw, unsigned win)
45 {
46  unsigned half_win=0.5*(win-1);
47  PolarScan<T> scan(raw);
48  Statistic<T> fit;
49  int pre,post;
50  for(unsigned i=0;i<raw.rows();i++)
51  {
52  fit.clear();
53  for(unsigned j=0;j<half_win+1;j++)
54  if(good(raw,i,j))
55  fit.feed(0*raw.cell_size,raw(i,j));
56  scan.set(i,0,fit.compute_slope());
57  if(scan(i,0)!=scan(i,0)) scan.set(i,0,raw.nodata); // changed to undetect
58 
59  for(unsigned j=1;j<raw.beam_size;j++)
60  {
61  pre=j-half_win-1;
62  post=j+half_win;
63  if(pre>=0)
64  if(good(raw,i,pre))
65  fit.slim(pre*raw.cell_size,raw(i,pre));
66  if(post<raw.beam_size)
67  if(good(raw,i,post))
68  fit.feed(post*raw.cell_size,raw(i,post));
69  scan.set(i,j,fit.compute_slope());
70  if(scan(i,j)!=scan(i,j)) scan.set(i,j,raw.nodata); // changed to undetect
71  }
72  }
73  return scan;
74 }
75 
83 template<typename T>
84 PolarScan<T> make_rms_scan(const PolarScan<T>& raw, unsigned len, unsigned wid)
85 {
86  unsigned half_len=0.5*(len-1);
87  unsigned half_wid=0.5*(wid-1);
88  PolarScan<T> scan(raw);
89  Statistic<T> rms;
90  int pre_l,post_l;
91  int pre_w,post_w;
92  for(unsigned i=0;i<raw.rows();i++)
93  {
94  rms.clear();
95  for(unsigned j=0;j<half_len+1;j++)
96  {
97  if(good(raw,i,j)) rms.feed(raw(i,j));
98  for(unsigned k=1;k<=half_wid;k++) // 0 è già inclusa, se half_wid=0 non dovrebbe partire il for
99  {
100  pre_w=i-k;
101  post_w=i+k;
102  if(pre_w<0)
103  pre_w=raw.beam_count+pre_w;
104  if(post_w>=raw.beam_count)
105  post_w=post_w-raw.beam_count;
106  if(good(raw,pre_w,j))
107  rms.feed(raw(pre_w,j));
108  if(good(raw,post_w,j))
109  rms.feed(raw(post_w,j));
110  }
111  }
112  scan.set(i,0,rms.compute_dev_std());
113  if(scan(i,0)!=scan(i,0)) scan.set(i,0,raw.nodata); // changed to undetect
114  for(unsigned j=1;j<raw.beam_size;j++)
115  {
116  pre_l=j-half_len-1;
117  post_l=j+half_len;
118  if(pre_l>=0)
119  {
120  if(good(raw,i,pre_l))
121  rms.slim(raw(i,pre_l));
122  for(unsigned k=1;k<=half_wid;k++)
123  {
124  pre_w=i-k;
125  post_w=i+k;
126  if(pre_w<0)
127  pre_w=raw.beam_count+pre_w;
128  if(post_w>=raw.beam_count)
129  post_w=post_w-raw.beam_count;
130  if(good(raw,pre_w,pre_l))
131  rms.slim(raw(pre_w,pre_l));
132  if(good(raw,post_w,pre_l))
133  rms.slim(raw(post_w,pre_l));
134  }
135  }
136  if(post_l<raw.beam_size)
137  {
138  if(good(raw,i,post_l))
139  rms.feed(raw(i,post_l));
140  for(unsigned k=1;k<=half_wid;k++)
141  {
142  pre_w=i-k;
143  post_w=i+k;
144  if(pre_w<0)
145  pre_w=raw.beam_count+pre_w;
146  if(post_w>=raw.beam_count)
147  post_w=post_w-raw.beam_count;
148  if(good(raw,pre_w,post_l))
149  rms.feed(raw(pre_w,post_l));
150  if(good(raw,post_w,post_l))
151  rms.feed(raw(post_w,post_l));
152  }
153  }
154  scan.set(i,j,rms.compute_dev_std());
155  if(scan(i,j)!=scan(i,j)) scan.set(i,j,raw.nodata);
156  }
157  }
158  return scan;
159 }
160 
168 template<typename T>
169 PolarScan<T> make_filter_scan(const PolarScan<T>& raw, unsigned len, unsigned wid)
170 {
171  unsigned half_len=0.5*(len-1);
172  unsigned half_wid=0.5*(wid-1);
173  PolarScan<T> scan(raw);
174  Statistic<T> filter;
175  int pre_l,post_l;
176  int pre_w,post_w;
177  for(unsigned i=0;i<raw.rows();i++)
178  {
179  filter.clear();
180  for(unsigned j=0;j<half_len+1;j++)
181  {
182  if(good(raw,i,j)) filter.feed(raw(i,j));
183  for(unsigned k=1;k<=half_wid;k++)
184  {
185  pre_w=i-k;
186  post_w=i+k;
187  if(pre_w<0)
188  pre_w=raw.beam_count+pre_w;
189  if(post_w>=raw.beam_count)
190  post_w=post_w-raw.beam_count;
191  if(good(raw,pre_w,j))
192  filter.feed(raw(pre_w,j));
193  if(good(raw,post_w,j))
194  filter.feed(raw(post_w,j));
195  }
196  }
197  scan.set(i,0,filter.compute_mean());
198  if(scan(i,0)!=scan(i,0)) scan.set(i,0,raw.nodata); // change to undetect
199  for(unsigned j=1;j<raw.beam_size;j++)
200  {
201  pre_l=j-half_len-1;
202  post_l=j+half_len;
203  if(pre_l>=0)
204  {
205  if(good(raw,i,pre_l))
206  filter.slim(raw(i,pre_l));
207  for(unsigned k=1;k<=half_wid;k++)
208  {
209  pre_w=i-k;
210  post_w=i+k;
211  if(pre_w<0)
212  pre_w=raw.beam_count+pre_w;
213  if(post_w>=raw.beam_count)
214  post_w=post_w-raw.beam_count;
215  if(good(raw,pre_w,pre_l))
216  filter.slim(raw(pre_w,pre_l));
217  if(good(raw,post_w,pre_l))
218  filter.slim(raw(post_w,pre_l));
219  }
220  }
221  if(post_l<raw.beam_size)
222  {
223  if(good(raw,i,post_l))
224  filter.feed(raw(i,post_l));
225  for(unsigned k=1;k<=half_wid;k++)
226  {
227  pre_w=i-k;
228  post_w=i+k;
229  if(pre_w<0)
230  pre_w=raw.beam_count+pre_w;
231  if(post_w>=raw.beam_count)
232  post_w=post_w-raw.beam_count;
233  if(good(raw,pre_w,post_l))
234  filter.feed(raw(pre_w,post_l));
235  if(good(raw,post_w,post_l))
236  filter.feed(raw(post_w,post_l));
237  }
238  }
239  scan.set(i,j,filter.compute_mean());
240  if(scan(i,j)!=scan(i,j)) scan.set(i,j,raw.nodata);
241  }
242  }
243  return scan;
244 }
245 
252 template<typename T>
253 PolarScan<T> make_gradient_azimuth_scan(const PolarScan<T>& raw)
254 {
255  PolarScan<T> scan(raw);
256 
257  for(unsigned rg=0;rg<raw.beam_size;rg++)
258  {
259  //cout<<"0"<<" "<<rg<<endl;
260  if(good(raw,raw.beam_count-1,rg) && good(raw,0,rg))
261  scan(0,rg) = raw(raw.beam_count-1,rg)-raw(0,rg);
262  else scan(0,rg) = 0.;
263  }
264  for(unsigned az=1;az<raw.beam_count;az++)
265  for(unsigned rg=0;rg<raw.beam_size;rg++)
266  {
267  //cout<<az<<" "<<rg<<endl;
268  if(good(raw,az-1,rg) && good(raw,az,rg) )
269  scan(az,rg) = raw(az-1,rg)-raw(az,rg);
270  else scan(az,rg) = 0.;
271  }
272  scan*=(0.5*raw.beam_count/M_PI);
273  return scan;
274 }
275 
283 template<typename T>
284 PolarScan<T> make_gradient_elevation_scan(const PolarScan<T>& low, const PolarScan<T>& up)
285 {
286  PolarScan<T> scan(up);
287  for(unsigned az=0;az<up.beam_count;az++)
288  for(unsigned rg=0;rg<up.beam_size;rg++)
289  {
290  if(good(up,az,rg) && good(low,az,rg))
291  scan(az,rg) = up(az,rg)-low(az,rg);
292  else scan(az,rg) = 0.;
293  }
294  scan*=abs(up.elevation-low.elevation)*M_PI/180.;
295  return scan;
296 }
297 
298 //===========================Volume manipulators=============================//
299 
300 namespace volume {
301 
309 template<typename T>
310 void moving_average_slope(const Volume<T>& raw, Volume<T>& vol, double slope_range, bool force_check_undetect=false)
311 {
312  unsigned window_size;
313  vol.clear();
314  check_undetect=force_check_undetect;
315  //this->quantity=raw.quantity.quantity_slope(); // TODO: è complesso ma si potrebbe
316  for(unsigned i=0;i<raw.size();i++)
317  {
318  window_size=1+2*std::floor(0.5*slope_range/raw[i].cell_size);
319  vol.push_back(make_slope_scan(raw[i],window_size));
320  }
321 }
322 
330 template<typename T>
331 void textureSD(const Scans<T>& raw, Scans<T>& vol, double filter_range, bool force_check_undetect=false)
332 {
333  textureSD(raw,vol,filter_range,0.,force_check_undetect);
334 }
335 
344 template<typename T>
345 void textureSD(const Scans<T>& raw, Scans<T>& vol, double filter_range, double filter_azimuth=0. ,bool force_check_undetect=false)
346 {
347  unsigned window_length;
348  unsigned window_width;
349  vol.clear();
350  check_undetect=force_check_undetect;
351  vol.quantity=raw.quantity;
352  vol.units=raw.units;
353  for(unsigned i=0;i<raw.size();i++)
354  {
355  window_length=1+2*std::floor(0.5*filter_range/raw[i].cell_size);
356  window_width=1+2*std::floor(0.5*filter_azimuth/(360./raw[i].beam_count));
357  vol.push_back(make_rms_scan(raw[i], window_length, window_width));
358  }
359 }
360 
368 template<typename T>
369 void filter(const Volume<T>& raw, Volume<T>& vol, double filter_range, bool force_check_undetect=false)
370 {
371  filter(raw,vol,filter_range,0.,force_check_undetect);
372 }
373 
382 template<typename T>
383 void filter(const Volume<T>& raw, Volume<T>& vol, double filter_range, double filter_azimuth=0., bool force_check_undetect=false)
384 {
385  unsigned window_length;
386  unsigned window_width;
387  vol.clear();
388  check_undetect=force_check_undetect;
389  vol.quantity=raw.quantity;
390  vol.units=raw.units;
391  for(unsigned i=0;i<raw.size();i++)
392  {
393  window_length=1+2*std::floor(0.5*filter_range/raw[i].cell_size);
394  window_width=1+2*std::floor(0.5*filter_azimuth/(360./raw[i].beam_count));
395  vol.push_back(make_filter_scan(raw[i], window_length, window_width));
396  }
397 }
398 
404 template<typename T>
405 void gradient_azimuth(const Volume<T>& raw, Volume<T>& vol, bool force_check_undetect=false)
406 {
407  vol.clear();
408  check_undetect=force_check_undetect;
409  vol.quantity=raw.quantity;
410  vol.units=raw.units;
411  for(unsigned el=0;el<raw.size();el++)
412  vol.push_back(make_gradient_azimuth_scan(raw[el]));
413 }
414 
420 template<typename T>
421 void gradient_elevation(const Volume<T>& raw, Volume<T>& vol, bool force_check_undetect=false)
422 {
423  vol.clear();
424  check_undetect=force_check_undetect;
425  vol.quantity=raw.quantity;
426  vol.units=raw.units;
427  vol.push_back(make_gradient_elevation_scan(raw[1],raw[0]));
428  for(unsigned el=1;el<raw.size();el++)
429  vol.push_back(make_gradient_elevation_scan(raw[el-1],raw[el]));
430 }
431 
437 template<typename T>
438 void lin2dB(Volume<T>& lin, Volume<T>& dB)
439 {
440  //this->quantity=lin.quantity.lin2dB(); // TODO: not yet implemented
441  dB.clear();
442  for(unsigned i=0;i<lin.size();i++)
443  {
444  dB.push_back(PolarScan<T>(lin[i].beam_count,lin[i].beam_size,0.));
445  dB[i].cell_size = lin[i].cell_size;
446  dB[i].elevation = lin[i].elevation;
447  dB[i].block(0,0,lin[i].beam_count,lin[i].beam_size)=lin[i].log10();
448  dB[i].array()*=10.;
449  }
450 }
451 
457 template<typename T>
458 void dB2lin(Volume<T>& dB, Volume<T>& lin)
459 {
460  //this->quantity=DB.quantity.dB2lin(); // TODO: not yet implemented
461  lin.clear();
462  for(unsigned i=0;i<dB.size();i++)
463  {
464  lin.push_back(PolarScan<T>(dB[i].beam_count,dB[i].beam_size,0.));
465  lin[i].cell_size = dB[i].cell_size;
466  lin[i].elevation = dB[i].elevation;
467  lin[i].block(0,0,dB[i].beam_count,dB[i].beam_size) = dB[i]*0.1;
468  lin[i].block(0,0,dB[i].beam_count,dB[i].beam_size) = lin[i].exp10();
469  }
470 }
471 
476 template<typename T>
477 void lin2dB(Volume<T>& lin)
478 {
479  //this->quantity=lin.quantity.lin2dB(); // TODO: not yet implemented
480  for(unsigned i=0;i<lin.size();i++)
481  lin[i].block(0,0,lin[i].beam_count,lin[i].beam_size)=lin[i].log10();
482  lin*=10.;
483 }
484 
489 template<typename T>
490 void dB2lin(Volume<T>& dB)
491 {
492  //this->quantity=DB.quantity.dB2lin(); // TODO: not yet implemented
493  dB*=0.1;
494  for(unsigned i=0;i<dB.size();i++)
495  dB[i].block(0,0,dB[i].beam_count,dB[i].beam_size) = dB.scan(i).exp10();
496 }
497 
498 } // namespace volume
499 } // namespace radarelab
Definisce le principali strutture che contengono i dati.
Class to manage statistical information about streamed data.