Ice Sheet System Model  4.18
Code documentation
Data Structures | Public Member Functions | Data Fields | Private Attributes
Quadtree Class Reference

#include <Quadtree.h>

Data Structures

class  QuadtreeBox
 

Public Member Functions

 Quadtree ()
 
 Quadtree (double xmin, double xmax, double ymin, double ymax, int maxdepth_in)
 
 ~Quadtree ()
 
void Add (Observation *observation)
 
void AddAndAverage (double x, double y, double value)
 
void ClosestObs (int *pindex, double x, double y)
 
void DeepEcho (void)
 
void Echo (void)
 
void IntergerCoordinates (int *xi, int *yi, double x, double y)
 
QuadtreeBoxNewQuadtreeBox (double xcenter, double ycenter, double length)
 
QuadtreeBoxNewQuadtreeBox (QuadtreeBox *master, int index)
 
void QuadtreeDepth (int *A, int xi, int yi)
 
void QuadtreeDepth2 (int *A, int xi, int yi)
 
void RangeSearch (int **pindices, int *pnobs, double x, double y, double range)
 

Data Fields

int MaxDepth
 
QuadtreeBoxroot
 
int NbQuadtreeBox
 
int NbObs
 

Private Attributes

DataSetboxcontainer
 

Detailed Description

Definition at line 7 of file Quadtree.h.

Constructor & Destructor Documentation

◆ Quadtree() [1/2]

Quadtree::Quadtree ( )

Definition at line 83 of file Quadtree.cpp.

83  {/*{{{*/
84  _error_("Constructor not supported");
85 
86 }

◆ Quadtree() [2/2]

Quadtree::Quadtree ( double  xmin,
double  xmax,
double  ymin,
double  ymax,
int  maxdepth_in 
)

Definition at line 88 of file Quadtree.cpp.

88  {/*{{{*/
89 
90  /*Intermediaries*/
91  double length;
92 
93  /*Initialize fields*/
94  this->MaxDepth=maxdepth;
95  this->NbQuadtreeBox=0;
96  this->NbObs=0;
97 
98  /*Create container*/
99  this->boxcontainer=new DataSet();
100 
101  /*Create Root, pointer toward the main box*/
102  length=max(xmax-xmin,ymax-ymin);
103  this->root=NewQuadtreeBox(xmin+length/2,ymin+length/2,length);
104 }

◆ ~Quadtree()

Quadtree::~Quadtree ( )

Definition at line 106 of file Quadtree.cpp.

106  {/*{{{*/
107 
108  delete boxcontainer;
109  root=NULL;
110 
111  }

Member Function Documentation

◆ Add()

void Quadtree::Add ( Observation observation)

Definition at line 115 of file Quadtree.cpp.

115  {/*{{{*/
116 
117  /*Intermediaries*/
118  int xi,yi,ij,level,levelbin;
119  QuadtreeBox **pbox = NULL; // pointer toward current box b
120  QuadtreeBox **pmaster = NULL; // pointer toward master of b
121  QuadtreeBox *box = NULL; // current box b
122  QuadtreeBox *slave = NULL; // suslaveox of b (if necessary)
123  Observation *obs[4];
124 
125  /*Get integer coodinates*/
126  xi = observation->xi;
127  yi = observation->yi;
128 
129  /*Initialize levels*/
130  level = 0;
131  levelbin = (1L<<this->MaxDepth);// = 2^30
132 
133  /*Get inital box (the largest)*/
134  pmaster = &root;
135  pbox = &root;
136 
137  /*Find the smallest box where the observation is located*/
138  while((box=*pbox) && (box->nbitems<0)){
139 
140  /*Go down one level (levelbin = 00100 -> 00010)*/
141  levelbin>>=1; level+=1; _assert_(level<this->MaxDepth);
142 
143  /*Get next box according to the bit value (levelbin)*/
144  pmaster = pbox;
145  pbox = &box->box[IJ(xi,yi,levelbin)];
146  }
147  _assert_(levelbin>0);
148 
149  /*Now, try to add the vertex, if the box is full (nbitems=4), we have to divide it in 4 new boxes*/
150  while((box=*pbox) && (box->nbitems==4)){
151 
152  /*Copy the 4 observation in the current Quadtreebox*/
153  obs[0] = box->obs[0];
154  obs[1] = box->obs[1];
155  obs[2] = box->obs[2];
156  obs[3] = box->obs[3];
157 
158  /*set nbitems as -1 (now holding boxes instead of observations)*/
159  box->nbitems = -1;
160  box->box[0] = NULL;
161  box->box[1] = NULL;
162  box->box[2] = NULL;
163  box->box[3] = NULL;
164 
165  /*Go down one level (levelbin = 00010 -> 00001)*/
166  levelbin>>=1; level+=1; _assert_(level<this->MaxDepth);
167 
168  /*Put the four observations in the new boxes*/
169  for (int k=0;k<4;k++){
170 
171  /*Get box for observation number k*/
172  ij = IJ(obs[k]->xi,obs[k]->yi,levelbin);
173  slave = box->box[ij];
174  if(!slave){
175  box->box[ij] = NewQuadtreeBox(box,ij);
176  slave = box->box[ij];
177  }
178  slave->obs[slave->nbitems++] = obs[k];
179  }
180 
181  /*Get the suslaveox where the current observation is located*/
182  ij = IJ(xi,yi,levelbin);
183  pmaster = pbox;
184  pbox = &box->box[ij];
185  }
186 
187  /*alloc the QuadtreeBox if necessary and add current observation*/
188  box = *pbox;
189  if(!box){
190  ij = IJ(xi,yi,levelbin);
191  box = *pbox = NewQuadtreeBox(*pmaster,ij);
192  }
193  box->obs[box->nbitems++]=observation;
194  NbObs++;
195 
196 }/*}}}*/

◆ AddAndAverage()

void Quadtree::AddAndAverage ( double  x,
double  y,
double  value 
)

Definition at line 197 of file Quadtree.cpp.

197  {/*{{{*/
198 
199  QuadtreeBox **pbox = NULL;
200  QuadtreeBox *box = NULL;
201  int xi,yi;
202  int levelbin;
203  int index;
204  double length,length2;
205 
206  /*Get integer coodinates*/
207  this->IntergerCoordinates(&xi,&yi,x,y);
208 
209  /*Initialize level*/
210  levelbin = (1L<<this->MaxDepth);// = 2^30
211 
212  /*Get inital box (the largest)*/
213  pbox=&root;
214 
215  /*Find the smallest box where this point is located*/
216  while((box=*pbox) && (box->nbitems<0)){
217  levelbin>>=1;
218  pbox = &box->box[IJ(xi,yi,levelbin)];
219  }
220 
221  /*Add obervation in this box (should be full)*/
222  if(box && box->nbitems==4){
223  index = 0;
224  length = pow(box->obs[0]->x - x,2) + pow(box->obs[0]->y - y,2);
225  for(int i=1;i<4;i++){
226  length2 = pow(box->obs[i]->x - x,2) + pow(box->obs[i]->y - y,2);
227  if(length2<length){
228  index = i;
229  length = length2;
230  }
231  }
232 
233  /*We found the closest observation, now average observation (do not change xi and yi to avoid round off errors*/
234  box->obs[index]->x = (box->obs[index]->weight*box->obs[index]->x + x)/(box->obs[index]->weight+1.);
235  box->obs[index]->y = (box->obs[index]->weight*box->obs[index]->y + y)/(box->obs[index]->weight+1.);
236  box->obs[index]->xi= int((box->obs[index]->weight*double(box->obs[index]->xi) + double(xi))/(box->obs[index]->weight+1.));
237  box->obs[index]->yi= int((box->obs[index]->weight*double(box->obs[index]->yi) + double(yi))/(box->obs[index]->weight+1.));
238  box->obs[index]->value = (box->obs[index]->weight*box->obs[index]->value + value)/(box->obs[index]->weight+1.);
239  box->obs[index]->weight += 1.;
240  }
241  else{
242  _error_("Box is not full");
243  }
244 }/*}}}*/

◆ ClosestObs()

void Quadtree::ClosestObs ( int *  pindex,
double  x,
double  y 
)

Definition at line 245 of file Quadtree.cpp.

245  {/*{{{*/
246 
247  QuadtreeBox **pbox = NULL;
248  QuadtreeBox *box = NULL;
249  int xi,yi;
250  int levelbin;
251  int index = -1;
252  double length,length2;
253 
254  /*Get integer coodinates*/
255  this->IntergerCoordinates(&xi,&yi,x,y);
256 
257  /*Initialize level*/
258  levelbin = (1L<<this->MaxDepth);// = 2^30
259 
260  /*Get inital box (the largest)*/
261  pbox=&root;
262 
263  /*Find the smallest box where this point is located*/
264  while((box=*pbox) && (box->nbitems<0)){
265  levelbin>>=1;
266  pbox = &box->box[IJ(xi,yi,levelbin)];
267  }
268 
269  /*Add obervation in this box (should be full)*/
270  if(box && box->nbitems>0){
271  index = box->obs[0]->index;
272  length = pow(box->obs[0]->x - x,2) + pow(box->obs[0]->y - y,2);
273  for(int i=1;i<box->nbitems;i++){
274  length2 = pow(box->obs[i]->x - x,2) + pow(box->obs[i]->y - y,2);
275  if(length2<length){
276  index = box->obs[i]->index;
277  length = length2;
278  }
279  }
280  }
281 
282  *pindex=index;
283 }/*}}}*/

◆ DeepEcho()

void Quadtree::DeepEcho ( void  )

Definition at line 284 of file Quadtree.cpp.

284  {/*{{{*/
285 
286  _printf_("Quadtree:\n");
287  _printf_(" MaxDepth = " << this->MaxDepth << "\n");
288  _printf_(" NbQuadtreeBox = " << this->NbQuadtreeBox << "\n");
289  _printf_(" NbObs = " << this->NbObs << "\n");
290  _printf_(" root = " << this->root << "\n");
291  boxcontainer->Echo();
292 
293 }/*}}}*/

◆ Echo()

void Quadtree::Echo ( void  )

Definition at line 294 of file Quadtree.cpp.

294  {/*{{{*/
295 
296  _printf_("Quadtree:\n");
297  _printf_(" MaxDepth = " << this->MaxDepth << "\n");
298  _printf_(" NbQuadtreeBox = " << this->NbQuadtreeBox << "\n");
299  _printf_(" NbObs = " << this->NbObs << "\n");
300  _printf_(" root = " << this->root << "\n");
301 
302 }/*}}}*/

◆ IntergerCoordinates()

void Quadtree::IntergerCoordinates ( int *  xi,
int *  yi,
double  x,
double  y 
)

Definition at line 303 of file Quadtree.cpp.

303  {/*{{{*/
304 
305  /*Intermediaries*/
306  double coefficient;
307  double xmin,ymin;
308 
309  /*Checks in debugging mode*/
310  _assert_(xi && yi);
311  _assert_(this->root);
312 
313  /*coeffIcoor is the coefficient used for integer coordinates:
314  * (x-xmin)
315  * xi = (2^30 -1) ---------
316  * length
317  * coefficient = (2^30 -1)/length
318  */
319  coefficient = double((1L<<this->MaxDepth)-1)/(this->root->length);
320  xmin = this->root->xcenter - this->root->length/2;
321  ymin = this->root->ycenter - this->root->length/2;
322 
323  *xi=int(coefficient*(x - xmin));
324  *yi=int(coefficient*(y - ymin));
325 }/*}}}*/

◆ NewQuadtreeBox() [1/2]

Quadtree::QuadtreeBox * Quadtree::NewQuadtreeBox ( double  xcenter,
double  ycenter,
double  length 
)

Definition at line 326 of file Quadtree.cpp.

326  {/*{{{*/
327 
328  /*Output*/
329  QuadtreeBox* newbox=NULL;
330 
331  /*Create and initialize a new box*/
332  newbox=new QuadtreeBox();
333  newbox->nbitems=0;
334  newbox->xcenter=xcenter;
335  newbox->ycenter=ycenter;
336  newbox->length=length;
337  newbox->box[0]=NULL;
338  newbox->box[1]=NULL;
339  newbox->box[2]=NULL;
340  newbox->box[3]=NULL;
341 
342  /*Add to container*/
343  this->boxcontainer->AddObject(newbox);
344  NbQuadtreeBox++;
345 
346  /*currentbox now points toward next quadtree box*/
347  return newbox;
348 }/*}}}*/

◆ NewQuadtreeBox() [2/2]

Quadtree::QuadtreeBox * Quadtree::NewQuadtreeBox ( QuadtreeBox master,
int  index 
)

Definition at line 349 of file Quadtree.cpp.

349  {/*{{{*/
350 
351  /*Output*/
352  QuadtreeBox* newbox=NULL;
353 
354  /*Checks in debugging mode*/
355  _assert_(master);
356 
357  /*Create and initialize a new box*/
358  newbox=new QuadtreeBox();
359  newbox->nbitems=0;
360  newbox->box[0]=NULL;
361  newbox->box[1]=NULL;
362  newbox->box[2]=NULL;
363  newbox->box[3]=NULL;
364  switch(index){
365  case 0:
366  newbox->xcenter=master->xcenter - master->length/4;
367  newbox->ycenter=master->ycenter - master->length/4;
368  break;
369  case 1:
370  newbox->xcenter=master->xcenter + master->length/4;
371  newbox->ycenter=master->ycenter - master->length/4;
372  break;
373  case 2:
374  newbox->xcenter=master->xcenter - master->length/4;
375  newbox->ycenter=master->ycenter + master->length/4;
376  break;
377  case 3:
378  newbox->xcenter=master->xcenter + master->length/4;
379  newbox->ycenter=master->ycenter + master->length/4;
380  break;
381  default:
382  _error_("Case " << index << " not supported");
383  }
384  newbox->length=master->length/2;
385 
386  /*Add to container*/
387  this->boxcontainer->AddObject(newbox);
388  NbQuadtreeBox++;
389 
390  /*currentbox now points toward next quadtree box*/
391  return newbox;
392 }/*}}}*/

◆ QuadtreeDepth()

void Quadtree::QuadtreeDepth ( int *  A,
int  xi,
int  yi 
)

Definition at line 410 of file Quadtree.cpp.

410  {/*{{{*/
411 
412  QuadtreeBox **pbox = NULL;
413  QuadtreeBox *box = NULL;
414  int level,levelbin;
415 
416  /*Initialize levels*/
417  level = 0;
418  levelbin = (1L<<this->MaxDepth);// = 2^30
419 
420  /*Get inital box (the largest)*/
421  pbox=&root;
422 
423  /*Find the smallest box where this point is located*/
424  while((box=*pbox) && (box->nbitems<0)){
425 
426  levelbin>>=1; level+=1; _assert_(level<this->MaxDepth);
427 
428  pbox = &box->box[IJ(xi,yi,levelbin)];
429  }
430  if(box && box->nbitems>0){
431  /*This box is not empty, add one level*/
432  level+=1;
433  }
434 
435  *A=level;
436 }/*}}}*/

◆ QuadtreeDepth2()

void Quadtree::QuadtreeDepth2 ( int *  A,
int  xi,
int  yi 
)

Definition at line 437 of file Quadtree.cpp.

437  {/*{{{*/
438 
439  QuadtreeBox **pbox = NULL;
440  QuadtreeBox *box = NULL;
441  int level,levelbin;
442 
443  /*Initialize levels*/
444  level = 0;
445  levelbin = (1L<<this->MaxDepth);// = 2^30
446 
447  /*Get inital box (the largest)*/
448  pbox=&root;
449 
450  /*Find the smallest box where this point is located*/
451  while((box=*pbox) && (box->nbitems<0)){
452 
453  levelbin>>=1; level+=1;
454 
455  pbox = &box->box[IJ(xi,yi,levelbin)];
456  }
457  if(box && box->nbitems>0){
458  /*This box is not empty, add one level*/
459  level+=1;
460  }
461 
462  /*If we were to add the vertex, get level*/
463  if(box && box->nbitems==4){
464  int ij;
465  bool flag=true;
466  while(flag){
467 
468  levelbin>>=1; level+=1;
469  if(level>this->MaxDepth){
470  level+=1;
471  break;
472  }
473 
474  /*loop over the four observations*/
475  ij=IJ(box->obs[0]->xi,box->obs[0]->yi,levelbin);
476  for (int k=1;k<4;k++){
477  if(IJ(box->obs[k]->xi,box->obs[k]->yi,levelbin) != ij){
478  flag = false;
479  }
480  }
481  if(IJ(xi,yi,levelbin)!=ij){
482  flag = false;
483  }
484  }
485  }
486 
487  *A=level;
488 }/*}}}*/

◆ RangeSearch()

void Quadtree::RangeSearch ( int **  pindices,
int *  pnobs,
double  x,
double  y,
double  range 
)

Definition at line 393 of file Quadtree.cpp.

393  {/*{{{*/
394 
395  /*Intermediaries*/
396  int nobs;
397  int *indices = NULL;
398 
399  /*Allocate indices (maximum by default*/
400  if(this->NbObs) indices = xNew<int>(this->NbObs);
401  nobs = 0;
402 
403  if(this->root) this->root->RangeSearch(indices,&nobs,x,y,range);
404 
405  /*Clean-up and return*/
406  *pnobs=nobs;
407  *pindices=indices;
408 
409 }/*}}}*/

Field Documentation

◆ boxcontainer

DataSet* Quadtree::boxcontainer
private

Definition at line 44 of file Quadtree.h.

◆ MaxDepth

int Quadtree::MaxDepth

Definition at line 47 of file Quadtree.h.

◆ root

QuadtreeBox* Quadtree::root

Definition at line 48 of file Quadtree.h.

◆ NbQuadtreeBox

int Quadtree::NbQuadtreeBox

Definition at line 49 of file Quadtree.h.

◆ NbObs

int Quadtree::NbObs

Definition at line 50 of file Quadtree.h.


The documentation for this class was generated from the following files:
Quadtree::boxcontainer
DataSet * boxcontainer
Definition: Quadtree.h:44
_assert_
#define _assert_(ignore)
Definition: exceptions.h:37
Observation
Definition: Observation.h:10
Quadtree::QuadtreeBox::ycenter
double ycenter
Definition: Quadtree.h:20
DataSet::AddObject
int AddObject(Object *object)
Definition: DataSet.cpp:252
_printf_
#define _printf_(StreamArgs)
Definition: Print.h:22
Quadtree::root
QuadtreeBox * root
Definition: Quadtree.h:48
Observation::yi
int yi
Definition: Observation.h:14
Quadtree::QuadtreeBox::obs
Observation * obs[4]
Definition: Quadtree.h:24
Observation::xi
int xi
Definition: Observation.h:14
DataSet::Echo
void Echo()
Definition: DataSet.cpp:307
Quadtree::MaxDepth
int MaxDepth
Definition: Quadtree.h:47
IJ
#define IJ(i, j, l)
Definition: Quadtree.cpp:79
Observation::index
int index
Definition: Observation.h:15
Quadtree::QuadtreeBox::length
double length
Definition: Quadtree.h:21
Quadtree::NewQuadtreeBox
QuadtreeBox * NewQuadtreeBox(double xcenter, double ycenter, double length)
Definition: Quadtree.cpp:326
Quadtree::QuadtreeBox::box
QuadtreeBox * box[4]
Definition: Quadtree.h:23
Quadtree::NbQuadtreeBox
int NbQuadtreeBox
Definition: Quadtree.h:49
Quadtree::QuadtreeBox::RangeSearch
void RangeSearch(int *indices, int *pnobs, double x, double y, double range)
Definition: Quadtree.cpp:534
_error_
#define _error_(StreamArgs)
Definition: exceptions.h:49
Quadtree::IntergerCoordinates
void IntergerCoordinates(int *xi, int *yi, double x, double y)
Definition: Quadtree.cpp:303
Quadtree::NbObs
int NbObs
Definition: Quadtree.h:50
max
IssmDouble max(IssmDouble a, IssmDouble b)
Definition: extrema.cpp:24
DataSet
Declaration of DataSet class.
Definition: DataSet.h:14
Quadtree::QuadtreeBox::xcenter
double xcenter
Definition: Quadtree.h:19