NASA Logo
Ocean Color Science Software

ocssw V2022
dim_mgr.cpp
Go to the documentation of this file.
1 /********************************************************************
2  dim_mgr.cpp - logic for the dimension manager
3  for a n dimensional array where you need only m (< n) dimensions
4  manage the storege for just the m-n dimensions for only the grid boxes
5  needed - this saves reading a big block in when you only need some
6  kind of slice.
7 
8  Modification history:
9  Programmer Date Description of change
10  ---------- ---- ---------------------
11  W. Robinson 8 Jan 2020 Original development
12  W. Robinson Nov 2020 make as a class in c++
13 
14 ********************************************************************/
15 #include "dim_mgr.hpp"
16 
17 // for class dim_mgr
18 // Separate the constructor out
19 
20  // ndim is # of dimensions to be managed
21 /********************************************************************
22  dim_mgr - 1st constructor form for dim_mgr class
23 
24  ndim I # of dimensions under management
25 ********************************************************************/
26  dim_mgr::dim_mgr(int32_t ndim) {
27  hash_tbl_siz = 7;
28  this->hash_tbl_siz = hash_tbl_siz;
29  this->const_hlp(ndim);
30  //cout << __FILE__ << ": 1-arg dim_mgr\n";
31  }
32 
33 /********************************************************************
34  dim_mgr - 2nd constructor form for dim_mgr class
35 
36  ndim I # of dimensions under management
37  hash_tbl_siz I # elements in the hash tablew
38 ********************************************************************/
39  dim_mgr::dim_mgr(int32_t ndim, int32_t hash_tbl_sz) {
40  hash_tbl_siz = hash_tbl_sz;
41  this->const_hlp(ndim);
42  this->hash_tbl_siz = hash_tbl_siz;
43  //cout << __FILE__ << ": 2-arg dim_mgr\n";
44  }
45 
46 /********************************************************************
47  const_hlp - helper for the 2 constructors
48  This does the real set-up of the dimension manager
49 
50  ndim I # of dimensions under management
51 ********************************************************************/
52  void dim_mgr::const_hlp(int32_t ndim) {
53  this->ndim = ndim;
54  prev_int = NULL;
55  n_tot_dat_blobs = 0;
56  dim_info = ( dim_info_struc ** )malloc
57  ( ndim * sizeof( dim_info_struc * ) );
58  //cout << __FILE__ << ": Setting the # dims\n";
59  // allocate all dim_info structs
60  for( int i = 0; i < ndim; i++ )
61  {
62  dim_info[i] = (dim_info_struc *) malloc( sizeof(dim_info_struc) );
63  dim_info[i]->dim_coords = NULL; // so we know if not filled
64  dim_info[i]->nvals = 0; // ditto
65  dim_info[i]->type = 0; // ditto
66  }
67  // to set the hash table
68  hash_tbl = new vector <hash_entry_struc>[hash_tbl_siz];
69  gpt_hash_tbl = new vector <gpt_hash_struc>[hash_tbl_siz];
70 
71  // to set up the point information structure
72  int ncorner = pow( 2, ndim );
73  pt_info.pt_status = (int32_t *)malloc( ncorner * sizeof(int32_t ) );
74  pt_info.pt_base_loc = (int32_t *)malloc( ndim * sizeof(int32_t ) );
75  pt_info.wt = (double *)malloc( ndim * sizeof(double) );
76  pt_info.wt_pt = (double *)malloc( ncorner * sizeof(double) );
77  pt_info.dat_ptrs = (void **)malloc( ncorner * sizeof(void *) );
78  }
79 
80 /********************************************************************
81  ~dim_mgr - destructor - needed for clean up of the class
82 ********************************************************************/
84  delete[] hash_tbl;
85  delete[] gpt_hash_tbl;
86  for( int i = 0; i < ndim; i++ )
87  {
88  if( dim_info[i]->type == 1 )
89  free( dim_info[i]->dim_coords );
90  free(dim_info[i]);
91  }
92  // dispense with point info
93  free( pt_info.pt_status );
94  free( pt_info.pt_base_loc );
95  free( pt_info.wt );
96  free( pt_info.wt_pt );
97  free( pt_info.dat_ptrs );
98  //
99  free( dim_info );
100  }
101 
102  // record an addition of points
103 /********************************************************************
104  add_pts - add to the total count of data blobs under management
105 
106  Returns:
107  void - none
108  Parameters: (in calling order)
109  Type Name I/O Description
110  ---- ---- --- -----------
111  int32_t nptadd I # grid points added
112 ********************************************************************/
113  void dim_mgr::add_pts( int32_t nptadd ) {
114  n_tot_dat_blobs += nptadd;
115  }
116 
117 /********************************************************************
118  get_tot_blobs - get total # data blobs being managed
119 
120  Returns:
121  int32_t # data blobs being managed
122 ********************************************************************/
123  int32_t dim_mgr::get_tot_blobs() { return n_tot_dat_blobs; }
124 
125 /********************************************************************
126  init_dim - initialize information for a dimension. First form is for
127  the dimension coordinates specified in an array
128 
129  Returns:
130  int32_t - status 0 = good
131 
132  Parameters: (in calling order)
133  Type Name I/O Description
134  ---- ---- --- -----------
135  int32_t dim_num I Dimension to set
136  int32_t nvals I # grid pts in this dim
137  double min I min of dimension range
138  double max I max of dimension range
139 ********************************************************************/
140  int32_t dim_mgr::init_dim( int32_t dim_num, int32_t nvals,
141  double min, double max ) {
142  // set the dim_info_struc if not set yet
143  if( dim_info == NULL ) {
144  cout << __FILE__ << __LINE__ << " dim_info array is null\n";
145  return 0;
146  }
147  if( dim_info[dim_num]->nvals != 0 ) {
148  cout << __FILE__ << __LINE__ << " dim_info array entry " << dim_num <<
149  " is filled already\n";
150  return 0;
151  }
152  dim_info[dim_num]->nvals = nvals;
153  dim_info[dim_num]->min = min;
154  dim_info[dim_num]->max = max;
155  dim_info[dim_num]->type = 0;
156  dim_info[dim_num]->delta = ( max - min ) / ( nvals - 1 );
157  return 0;
158  }
159 
160 /********************************************************************
161  init_dim - initialize information for a dimension. Second form is for
162  the dimension coordinates specified in an array
163 
164  Parameters: (in calling order)
165  Type Name I/O Description
166  ---- ---- --- -----------
167  int32_t dim_num I Dimension to set
168  int32_t nvals I # grid pts in this dim
169  double * dim_coords I set of coordinates for this
170  dimension
171 ********************************************************************/
172  int32_t dim_mgr::init_dim( int32_t dim_num, int32_t nvals,
173  double *dim_coords ) {
174 
175  if( dim_info == NULL ) {
176  cout << __FILE__ << __LINE__ << " dim_info array is null\n";
177  return 0;
178  }
179  if( dim_info[dim_num]->dim_coords != NULL ) {
180  cout << __FILE__ << __LINE__ << " dim_info array entry " << dim_num <<
181  " is filled already\n";
182  return 0;
183  }
184  dim_info[dim_num]->nvals = nvals;
185  dim_info[dim_num]->min = dim_coords[0];
186  dim_info[dim_num]->max = dim_coords[ nvals - 1 ];
187  dim_info[dim_num]->type = 1;
188  dim_info[dim_num]->dim_coords =
189  (double *)malloc( nvals * sizeof( double ) );
190  for( int i = 0; i < nvals; i++ )
191  dim_info[dim_num]->dim_coords[i] = dim_coords[i];
192  return 0;
193  }
194 
195 /********************************************************************
196  update_new_data - Update info about new data under management
197 
198  This will refresh the interval's data status and data pointer using
199  the pt_info structure the user updated.
200 
201 ********************************************************************/
203  int32_t npt = pow( 2, ndim );
204 
205  // just may as well update all points
206  // WDR 27 Aug 2024 do not update if prev_int is NULL
207  // to test use of prev_int
208  if( prev_int != (interval_struc *)NULL )
209  {
210  for( int i = 0; i < npt; i++ )
211  {
212  prev_int->dat_info[i]->dat_status = pt_info.pt_status[i];
213  prev_int->dat_info[i]->dat = pt_info.dat_ptrs[i];
214  }
215  }
216  }
217 
218 /********************************************************************
219  mng_pt - for a point in a grid, find or set up the location for the
220  associated data and return that and the weights to apply
221 
222  Returns:
223  pt_info_struc * the structure with all point information
224  needed to interpolagte the data
225 
226  Parameters: (in calling order)
227  Type Name I/O Description
228  ---- ---- --- -----------
229  double * pt I coordinates of the point to
230  locate
231  int32_t access_id I A last time of access to
232  assign to the bottom interval
233  the lower the value, the
234  more time since last access
235  int32_t * status O return status 0 = good,
236  1 = code failure, you
237  shoiuld fail at this point
238  2 = not valid point (not in
239  the grid bounds) treat as
240  unable to do the point
241  Modification history:
242  Programmer Date Description of change
243  ---------- ---- ---------------------
244  W. Robinson 1 Sep 2020 based loosly on sparse_mng_pt()
245 ********************************************************************/
246  pt_info_struc * dim_mgr::mng_pt( double *pt, int32_t access_id,
247  int32_t *status ) {
248 
249  int32_t idim, idat, ndat, same_interval;
250  int32_t *ix, *ix_real, *offset_lcl, acc_mode, ret;
251  interval_struc *found_int;
252 
253  *status = 0; // start with good status
254  /*
255  * If the new point is the same as the old, just return the same point info
256  */
257  // not in std c++ so... if( isnormal( old_pt[0] )
258  if( old_pt[0] == old_pt[0] )
259  {
260  int32_t match = 1;
261  for( idim = 0; idim < ndim; idim++ )
262  {
263  if( *( pt + idim ) != *( old_pt + idim ) )
264  {
265  match = 0;
266  break;
267  }
268  }
269  if( match == 1 )
270  {
271  *status = old_status;
272  return &pt_info;
273  }
274  }
275  for( idim = 0; idim < ndim; idim++ )
276  *( old_pt + idim ) = *( pt + idim );
277  /*
278  * set up some local storage
279  */
280  if( ( ( ix = (int32_t *) malloc( ndim * sizeof( int32_t ) ) ) == NULL ) ||
281  ( ( ix_real = (int32_t *)
282  malloc( ndim * sizeof( int32_t ) ) ) == NULL ) ||
283  ( ( offset_lcl = (int32_t *)
284  malloc( ndim * sizeof( int32_t ) ) ) == NULL ) )
285  {
286  printf( "%s, %d: Unable to allocate the dim index array\n", __FILE__,
287  __LINE__ );
288  *status = 1;
289  return &pt_info;
290  }
291  /*
292  * also set the weight and point weight if needed
293  */
294  ndat = pow( 2., ndim );
295  /*
296  * get the index of the point in all dims and the weight to it
297  */
298  if( sparse_get_loc( dim_info, ndim, pt, ix, pt_info.wt, pt_info.wt_pt )
299  != 0 )
300  {
301  *status = 2;
302  old_status = 2;
303  return &pt_info;
304  }
305  /*
306  printf( "\n\n%s - index: %d,%d, weight: %f,%f\n", __FILE__,
307  ix[0], ix[1], pt_info.wt[0], pt_info.wt[1] );
308  printf( "data point: (%f, %f)\n", pt[0], pt[1] );
309  */
310  /*
311  * Its likely that the previous point is in the same grid box as the current
312  * point. In this case, just return the previously found interval
313  */
314  same_interval = 0;
315  pt_info.interval_needs_data = 0;
316  if( prev_int != NULL )
317  {
318  /* see if we are in the same interval */
319  same_interval = 1;
320  for( idim = 0; idim < ndim; idim++ )
321  {
322  if( ix[idim] != prev_int->dat_info[0]->ix_arr[idim] )
323  {
324  same_interval = 0;
325  break;
326  }
327  }
328  }
329 
330  if( same_interval == 1 )
331  {
332  found_int = prev_int;
333  found_int->access_id = access_id;
334  }
335  else
336  {
337  /*
338  * get/set the interval for this grid point
339  */
340  acc_mode = 0;
341  valarray<int> ix_val(ndim);
342  for( int i = 0; i < ndim; i++ )
343  ix_val[i] = ix[i];
344  found_int = access( ix_val, acc_mode );
345 
346  /*
347  * Set the found interval's access id to the one passed in
348  */
349  found_int->access_id = access_id;
350  /*
351  * if interval is new, we need to donate any data areas shared
352  * to the new interval
353  */
354  if( found_int->dat_filled == 0 )
355  {
356  idim = -1; // set this way to start at the top
357 
358  if( share_gpt( found_int, ix ) != 0 ) {
359  *status = 1;
360  return &pt_info;
361  }
362  // this should be obsolete
363  /* set the interval filled value to yes */
364  found_int->dat_filled = 1;
365  /* loop through all data pointers in the new interval */
366  /* and set up the data structures, less the pointer to the data */
367  for( idat = 0; idat < ndat; idat++ )
368  {
369  /* if the data for this corner is not there, create and fill it */
370  if( found_int->dat_info[idat] == NULL )
371  {
372  valarray<int> ixv_off(ndim);
373  pt_info.interval_needs_data = 1; // flag that data needed
374  // for this interval / grid box
375  /* set up each data info struct */
376  found_int->dat_info[idat] =
377  (data_info_struc *) malloc( sizeof( data_info_struc ) );
378  /* get index for this point */
379  linear_to_offset( ndim, idat, offset_lcl );
380  found_int->dat_info[idat]->ix_arr = (int32_t *) malloc( ndim *
381  sizeof( int32_t ) );
382  for( idim = 0; idim < ndim; idim++ )
383  {
384  ix_real[idim] = ix[idim] + offset_lcl[idim];
385  ixv_off[idim] = ix_real[idim];
386  found_int->dat_info[idat]->ix_arr[idim] = ix_real[idim];
387  }
388  /*
389  * set data status to not filled and # intervals pointing to it to 1
390  */
391  found_int->dat_info[idat]->dat_status = -1; /* not filled here */
392  found_int->dat_info[idat]->n_intervals = 1;
393  /*
394  * add new data_info to the grid point hash table
395  */
396  ret = gpt_add( ixv_off, found_int->dat_info[idat] );
397  if( ret != 0 )
398  {
399  *status = 1;
400  return &pt_info;
401  }
402  }
403  else /* for data already there, incriment access count */
404  {
405  // WDR keep in case some functions needed here
406  }
407  }
408  }
409  /* remember the found interval and its dim location */
410  prev_int = found_int;
411  }
412 
413  /* end set up of (possible) new grid interval */
414 
415  free( ix ); free( ix_real ); free( offset_lcl );
416 
417  // put the needed info into the pt_info struct
418  for( int i = 0; i < ndim; i++ )
419  pt_info.pt_base_loc[i] = found_int->dat_info[0]->ix_arr[i];
420  for( int i = 0; i < ndat; i++ ) {
421  pt_info.pt_status[i] = found_int->dat_info[i]->dat_status;
422  // pt_info.wt_pt, and wt are already set
423  pt_info.dat_ptrs[i] = found_int->dat_info[i]->dat;
424  }
425  /*
426  * return success
427  */
428  *status = 0;
429  old_status = 0;
430  return &pt_info;
431  }
432  // all the rest of routines that go with mng_pt
433 
434  int32_t dim_mgr::sparse_get_loc( dim_info_struc **dim_info, int32_t ndim,
435  double *pt, int32_t *ix, double *wt, double *wt_pt )
436 /********************************************************************
437  sparse_get_loc - get the location of a point in all grid dimensions
438  and the weights
439 
440  Returns:
441  int32_t - status 0 = good
442 
443  Parameters: (in calling order)
444  Type Name I/O Description
445  ---- ---- --- -----------
446  dim_info_struc ** dim_info I array of descriptions of each
447  dimension
448  that is handled
449  int32_t ndim I # dimensions handled
450  double * pt I coordinates of the point to
451  locate
452  int32_t * ix O index of the grid interval
453  for all dimensions
454  double * wt O for each dimension, the weight
455  to apply to the index point
456  (1 - wt for the next index
457  point)
458  double * wt_pt O set of consolidated weights to
459  apply to the data blobs pointed
460  to by the dat_info_struc
461  Modification history:
462  Programmer Date Description of change
463  ---------- ---- ---------------------
464  W. Robinson 26 Nov 2019 Original development
465 
466 ********************************************************************/
467  {
468  int32_t idim, iint, full_len, msk;
469  /*
470  * we start with all composite weights at 1
471  */
472  full_len = pow( 2, ndim );
473  for( idim = 0; idim < full_len; idim++ )
474  *( wt_pt + idim ) = 1.;
475  /*
476  loop through all dimensions and get the point and weight
477  */
478  for( idim = 0; idim < ndim; idim++ )
479  {
480  /* for now, assume the grid is uneven and we look through each interval */
481  if( ( pt[idim] < dim_info[idim]->min ) ||
482  ( pt[idim] > dim_info[idim]->max ) )
483  {
484  //printf( "%s, %d - point is outside grid dimension range\n",
485  // __FILE__, __LINE__ );
486  return 1;
487  }
488  for( iint = 1; iint < dim_info[idim]->nvals; iint++ )
489  {
490  if( pt[idim] <= dim_info[idim]->dim_coords[iint] )
491  {
492  ix[idim] = iint - 1;
493 
494  wt[idim] = ( pt[idim] - dim_info[idim]->dim_coords[iint-1] ) /
495  ( dim_info[idim]->dim_coords[iint] -
496  dim_info[idim]->dim_coords[iint-1] );
497 
498  break;
499  }
500  }
501  /*
502  * This will modify the composite weight so it applies to each
503  * corner's data - apply the weight to the end point
504  * of that grid else 1 - weight
505  */
506  msk = pow( 2, idim );
507  for( iint = 0; iint < full_len; iint++ )
508  {
509  if( ( iint & msk ) == 0 ) /* apply 1 - wt */
510  *( wt_pt + iint ) *= 1. - wt[idim];
511  else
512  *( wt_pt + iint ) *= wt[idim];
513  }
514  }
515  return 0;
516  }
517 
518 interval_struc * dim_mgr::access( valarray<int>& s, int32_t acc_mode )
519 /********************************************************************
520  access search / insert interval defining the data at the grid point
521 
522  Returns:
523  interval_struc * pointer to the interval information
524 
525  Parameters: (in calling order)
526  Type Name I/O Description
527  ---- ---- --- -----------
528  valarray<int> s I set of coordinates of the point
529  int32_t acc_mode I action: 0 - find the interval,
530  if not there, make it. return
531  the interval pointer, 1 -
532  find interval and return NULL
533  if not found, 2 - remove an
534  entry based on coordinates s
535  Modification history:
536  Programmer Date Description of change
537  ---------- ---- ---------------------
538  W. Robinson 31 Aug 2020 Original development
539 ********************************************************************/
540  {
541  int32_t ndat = pow( 2, ndim );
542 
543  //Compute the index by using the hash function
544  int index = hash_func( s, hash_tbl_siz );
545  //Search the linked list at that specific index for the index array
546  for( long unsigned int i = 0; i < hash_tbl[index].size(); i++)
547  {
548  // the == for the 2 valarrays make a valarray of booleans,
549  // basically all dims match (all 1), so min is not = 0
550  // if all dimensions match...
551  if( ( hash_tbl[index][i].ix_arr == s ).min() != 0 )
552  {
553  // interval is found, we can either return interval info or erase it
554  if( acc_mode == 2 ) // remove NOTE that underlying structures
555  // not handled here!!! so fix before use
556  {
557  hash_tbl[index].erase( hash_tbl[index].begin() + i );
558  cout << "Set is removed" << endl;
559  return (interval_struc *) NULL;
560  }
561  else // return interval struc for the interval
562  {
563  //cout << "Set is found!" << endl;
564  //cout << "position: " << i << endl;
565  return hash_tbl[index][i].int_str;
566  }
567  }
568  }
569  // if not found,
570  //cout << "Set is not found!" << endl;
571  if( acc_mode == 0 ) // add a new entry
572  {
573  // create the new entry
574  hash_entry_struc hash_entry;
576 
577  insert = (interval_struc *)malloc( sizeof(interval_struc) );
578  insert->dat_filled = 0;
579  insert->access_id = 0; // tenp setting, done later
580  insert->dat_info = (data_info_struc **)
581  calloc( ndat, sizeof(data_info_struc *) );
582  hash_entry.ix_arr = s;
583  hash_entry.int_str = insert;
584  // Insert the element in the linked list at the particular hash index
585 
586  hash_tbl[index].push_back(hash_entry);
587  //cout << "Inserting new interval, hash item\n";
588  return insert;
589  }
590  else
591  {
592  return (interval_struc *) NULL;
593  }
594  cout << "Code should not get here!/n";
595  return (interval_struc *)NULL;
596  }
597 
598 // New - add a data info structure to the grid point hash table
599 
600 int32_t dim_mgr::gpt_add( valarray<int> s, data_info_struc *dat_info )
601 /*-------------------------------------------------------------------
602  gpt_add - add a data info structure to the grid point hash table
603 
604  Returns: int 0 all OK, 1 trying to insert to a existing coordinate
605 
606  Parameters: (in calling order)
607  Type Name I/O Description
608  ---- ---- --- -----------
609  valarray<int> s I set of coordinates of the point
610  data_info_struc * dat_info I descriptor for the data at the
611  grid point
612 
613  Modification history:
614  Programmer Date Description of change
615  ---------- ---- ---------------------
616  W. Robinson 30 Nov 2021 replace interval searching with this
617 
618 -------------------------------------------------------------------*/
619  {
620  // as with access, Compute the index by using the hash function
621  // and check for this coordinate already (there should not be one)
622  int index = hash_func( s, hash_tbl_siz );
623 
624  for( long unsigned int i = 0; i < gpt_hash_tbl[index].size(); i++)
625  {
626  if( ( gpt_hash_tbl[index][i].ix_arr == s ).min() != 0 )
627  {
628  // we have a match, which should not happen as this is insert
629  cout << __FILE__ << __LINE__ <<
630  "Error: grid point is already in hash table, Exiting" << endl;
631  return 1;
632  }
633  }
634  // Normal branch, insert the new hash entry
635  gpt_hash_struc hash_entry;
636 
637  hash_entry.ix_arr = s;
638  hash_entry.dat_info = dat_info;
639 
640  gpt_hash_tbl[index].push_back( hash_entry );
641  /*
642  cout << __FILE__ << ", " << __LINE__ << ", " <<
643  "Just added a grid point to hash" << endl;
644  */
645  return 0;
646  }
647 
648 int dim_mgr::hash_func(valarray<int> s, int32_t hash_tbl_siz )
649 /********************************************************************
650  hash_func - from the grid coordinates, make a hash entry number
651  This is one of the more odd parts of using a hashtable - get a good
652  distribution of table entries or the hash is not as useful
653 
654  Returns:
655  interval_struc * pointer to the interval information
656 
657  Parameters: (in calling order)
658  Type Name I/O Description
659  ---- ---- --- -----------
660  valarray<int> s I set of coordinates of the point
661  int32_t hash_tbl_siz I size of hash table
662 
663  Modification history:
664  Programmer Date Description of change
665  ---------- ---- ---------------------
666  W. Robinson 31 Aug 2020 Original development
667 ********************************************************************/
668  {
669  int retval;
670  long sum = 0;
671 
672  if( hash_ifirst == 0 ) {
673  int32_t ndim = s.size();
674  hash_ifirst = 1;
675  hash_h_mult = hash_tbl_siz / ( 2 * ndim );
676  hash_h_mult = hash_h_mult - ndim / 2;
677  }
678 
679  for( long unsigned int i = 0; i < s.size(); i++ )
680  sum += ( hash_h_mult + i ) * s[i];
681  sum = sum % hash_tbl_siz;
682 
683  retval = sum;
684  return retval;
685  }
686 
687 int32_t dim_mgr::share_gpt( interval_struc *intvl, int32_t *ix )
688 /********************************************************************
689  share_gpt will find grid points that are already set up and share
690  them with the current interval
691 
692  Returns:
693  int32_t - 0 if all good
694 
695  Parameters: (in calling order)
696  Type Name I/O Description
697  ---- ---- --- -----------
698  interval_struc * intvl I/O a new interval structure to
699  find data for
700  int32_t * ix I start grid location of the
701  interval
702  Modification history:
703  Programmer Date Description of change
704  ---------- ---- ---------------------
705  W. Robinson 1 Dec 2021 to hopefully simplify the grid point
706  data sharing process, replacing the
707  iterative int_share_data
708 ********************************************************************/
709  {
710  int32_t *ix_off, i;
711  valarray<int> s(ndim);
712  int32_t npt = pow( 2, ndim );
713  ix_off = ( int32_t * ) malloc( ndim * sizeof(int32_t) );
714  //
715  // look at all the corners of this interval
716  for( int ipt = 0; ipt < npt; ipt++ )
717  {
718  // get the grid point location from base and offset to this point
719  linear_to_offset( ndim, ipt, ix_off );
720  for( i = 0; i < ndim; i++ )
721  {
722  ix_off[i] += ix[i];
723  s[i] = ix_off[i];
724  }
725 
726  // see if the point is in the grid point hash already
727  int index = hash_func( s, hash_tbl_siz );
728  for( i = 0; i < (int32_t)gpt_hash_tbl[index].size(); i++ )
729  {
730  if( ( gpt_hash_tbl[index][i].ix_arr == s ).min() != 0 )
731  {
732  // found matching pt it in hash table, point to it from the interval
733  intvl->dat_filled = 1;
734  intvl->dat_info[ipt] = gpt_hash_tbl[index][i].dat_info;
735  intvl->dat_info[ipt]->n_intervals++; // one more ref to the point
736  break;
737  }
738  }
739  // new points are handled in mng_pt
740  }
741  free( ix_off );
742  return 0;
743  }
744 
745 /********************************************************************
746  purge will completely remove all entries from the dimension manager
747 
748  Returns:
749  int32_t - 0 all good
750 
751 ********************************************************************/
752 int32_t dim_mgr::purge()
753  {
754  int32_t nhash = hash_tbl_siz;
755  int32_t ndat = pow( 2, ndim );
756  // go to hash entry
757  // loop over each hash bin
758  for( int ihash = 0; ihash < nhash; ihash++ )
759  {
760  int32_t nentry = hash_tbl[ihash].size();
761  for( int iint = 0; iint < nentry; iint++ )
762  {
763  // go to interval struc
764  hash_entry_struc hash_ent = hash_tbl[ihash][0];
765 
766  interval_struc *int_str = hash_ent.int_str;
767 
768  // go to data structures
769  for( int idat = 0; idat < ndat; idat++ )
770  {
771  // rid data if # intervals accessing is 1
772  if( int_str->dat_info[idat]->n_intervals > 1 )
773  {
774  int_str->dat_info[idat]->n_intervals--;
775  }
776  else
777  {
778  // no more intervals point to this data struct,
779  // so free the data info structure
780  free( int_str->dat_info[idat]->dat );
781  free( int_str->dat_info[idat]->ix_arr );
782  free( int_str->dat_info[idat] );
783  }
784  }
785 
786  // free the array with data info structure addresses
787  free( int_str->dat_info );
788 
789  // free the interval struc
790  free( int_str );
791  // rid the top hash entry for bin i
792  hash_tbl[ihash].erase(hash_tbl[ihash].begin() );
793  }
794  }
795 
796  // Clean up all the grid point hash table entries
797  // the dat_info is already gone from the interval purge above, so
798  // the only thing left is all the entries in each gpt_hash
799  for( int ihash = 0; ihash < nhash; ihash++ )
800  {
801  int32_t nentry = gpt_hash_tbl[ihash].size();
802  for( int iint = 0; iint < nentry; iint++ )
803  {
804  gpt_hash_tbl[ihash].erase( gpt_hash_tbl[ihash].begin() );
805  }
806  }
807 
808  // the previous interval pointer is invalid, say so
809  prev_int = NULL;
810  n_tot_dat_blobs = 0;
811  // The rest can stay - the destructor handles that (I think)
812  return 0;
813  }
814 
815 /********************************************************************
816  prune - Remove selected intervals from management
817 
818  Returns: int32_t - 0 all good
819  access_id I erase intervals with access_id lower than access_id
820 ********************************************************************/
821 int32_t dim_mgr::prune(int32_t access_id ) {
822 
823  // start at dim_mgr
824  int32_t nhash = hash_tbl_siz;
825  int32_t ndat = pow( 2, ndim );
826 
827  // go to hash entry, loop over each hash bin
828  for( int ihash = 0; ihash < nhash; ihash++ )
829  {
830  // prune from list end to start so next entry location is unchanged
831  // even if an entry is deleted
832  int32_t nentry = hash_tbl[ihash].size();
833  for( int iint = ( nentry - 1 ); iint >= 0; iint-- )
834  {
835  // go to interval struc
836  hash_entry_struc hash_ent = hash_tbl[ihash][iint];
837 
838  interval_struc *int_str = hash_ent.int_str;
839  if( int_str->access_id < access_id )
840  {
841  // go to data structures
842  for( int idat = 0; idat < ndat; idat++ )
843  {
844  // rid data if # intervals accessed is 1 else lower interval count
845  if( int_str->dat_info[idat]->n_intervals > 0 )
846  {
847  int_str->dat_info[idat]->n_intervals--;
848  }
849  else
850  {
851  // no more intervals point to this data struct,
852  // leave the data info struc around to be handled when
853  // cleaning up the grid point hash
854  }
855  }
856  // WDR *** MAY need to free( int_str->dat_info );
857  free( int_str->dat_info );
858  // free the interval struc
859  free( int_str );
860  // rid the hash entry for hash ihash, list entry iint
861  hash_tbl[ihash].erase(hash_tbl[ihash].begin() + iint);
862  }
863  }
864  }
865  // go through the grid point hash table and remove entries with
866  // no accesses
867  for( int ihash = 0; ihash < nhash; ihash++ )
868  {
869  // prune from list end to start as before
870  int32_t nentry = gpt_hash_tbl[ihash].size();
871  for( int iint = ( nentry - 1 ); iint >= 0; iint-- )
872  {
873  // go to data info struc
874  gpt_hash_struc hash_ent = gpt_hash_tbl[ihash][iint];
875  data_info_struc *lcl_info = hash_ent.dat_info;
876  if( lcl_info->n_intervals == 0 )
877  {
878  // as no intervals point to this data, remove it
879  free( lcl_info->dat );
880  free( lcl_info->ix_arr );
881  free( lcl_info );
882  n_tot_dat_blobs--; // one less managed data element
883  gpt_hash_tbl[ihash].erase( gpt_hash_tbl[ihash].begin() + iint );
884  }
885  }
886  }
887  // This is much like purge, but should work a bit slower
888  prev_int = NULL; // resetprevious interval in case we removed it
889  return 0;
890  }
891 
892 /********************************************************************
893  dump_mgr - report all info on all managed intervals
894 
895  double pt I depth of reporting: 0 do all, 1 - leave out interval info
896 ********************************************************************/
897 void dim_mgr::dump_mgr( double *pt ) {
898  // print # bins in hash table
899  int32_t idepth = (int32_t)pt[0];
900  printf( "\nCurrent manager state:\n" );
901  printf( " # of hash tbl bins: %d, # data blobs managed: %d\n",
902  hash_tbl_siz, n_tot_dat_blobs );
903 
904  // loop all hash bins
905  for( int ibin = 0; ibin < hash_tbl_siz; ibin++ ) {
906  // print bin # an entries inbin
907  int32_t n_entry = hash_tbl[ibin].size();
908  printf( " bin # %d, # entries: %d\n", ibin, n_entry );
909  // loop bin entries
910  for( int ient = 0; ient < n_entry; ient++ ) {
911  // print entry # and interval index associated with this entry
912  printf( " Entry %d, interval index for this entry:\n ",
913  ient );
914  for( int ix = 0; ix < ndim; ix++ )
915  printf( " %d,", hash_tbl[ibin][ient].ix_arr[ix] );
916  printf( "\n" );
917  if( idepth < 1 ) {
918  printf( "-------------------------\n" );
919  printf( "Interval Information:\n" );
920  // call the interval dumper
921  dump_interval( hash_tbl[ibin][ient].int_str );
922  printf( "-------------------------\n" );
923  }
924  }
925  }
926  // Add a report for the grid point hash bins
927  printf( "\n\nGrid point hash table summary:\n" );
928  for( int ibin = 0; ibin < hash_tbl_siz; ibin++ )
929  {
930  // print bin # an entries inbin
931  int32_t n_entry = gpt_hash_tbl[ibin].size();
932  printf( "-------------------------\n" );
933  printf( " bin # %d, # entries: %d\n", ibin, n_entry );
934  // loop bin entries
935  for( int ient = 0; ient < n_entry; ient++ ) {
936  // print entry #
937  printf( " Entry #: %d, Point coords follow\n", ient );
938  for( int ix = 0; ix < ndim; ix++ )
939  printf( " %d,", gpt_hash_tbl[ibin][ient].ix_arr[ix] );
940  printf( "\n" );
941  printf( "data info Point coords\n" );
942  for( int ix = 0; ix < ndim; ix++ )
943  printf( " %d,", gpt_hash_tbl[ibin][ient].dat_info->ix_arr[ix] );
944  printf( "\n" );
945  printf( "data info # intervals pointed to: %d\n",
946  gpt_hash_tbl[ibin][ient].dat_info->n_intervals );
947  printf( "data info status: %d\n",
948  gpt_hash_tbl[ibin][ient].dat_info->dat_status );
949  printf( "data info address: %ld\n",
950  (long int)gpt_hash_tbl[ibin][ient].dat_info->dat );
951  }
952  }
953  }
954 
955 /********************************************************************
956  dump_last_interval - report information about last interval
957  = a grid box's info
958 
959  interval_struc * interval I Interval structure for the
960 
961 ********************************************************************/
962 void dim_mgr::dump_last_interval() { dump_interval( prev_int ); }
963 
964 /********************************************************************
965  dump_interval - report information about an interval = a grid box's info
966 
967  interval_struc * interval I Interval structure for the
968 
969 ********************************************************************/
970 void dim_mgr::dump_interval( interval_struc *interval ) {
971  if( interval != NULL ) {
972  int32_t npt = pow( 2, ndim );
973  printf( "Interval address: %ld access ID: %d\n",(long int)interval,
974  interval->access_id );
975  printf( " Grid loc (%d dims): ", ndim );
976  for( int i = 0; i < ndim; i++ )
977  printf( " %d", interval->dat_info[0]->ix_arr[i] );
978  printf( "\n" );
979  for( int i = 0; i < npt; i++ )
980  {
981  printf(
982  "# %d data blob, address: %ld, # intervals pointing to it: %d\n", i,
983  (long int)interval->dat_info[i]->dat,
984  interval->dat_info[i]->n_intervals );
985  printf( " status: %d, grid point indexes: \n ",
986  interval->dat_info[i]->dat_status );
987  for( int j = 0; j < ndim; j++ )
988  printf( "%d, ", interval->dat_info[i]->ix_arr[j] );
989  printf( "\n" );
990  }
991  } else {
992  printf( "Sorry, last interval is NULL (no intervals added yet, " );
993  printf( "or purge or prune done)\n" );
994  }
995  }
996 
997  // information reporting
998  int dim_mgr::get_ndim() { return ndim; }
999  int dim_mgr::get_hdim() { return hash_tbl_siz; }
1000 
1001 // The linear_to_offset does not need to be in the class, in fact, it
1002 // himders use from outside. So make it a non-member function here
1003 
1004 /********************************************************************
1005  linear_to_offset will convert a linear location in the data pointer
1006  array to a offset array
1007 
1008  Returns: int32_t - 0 if all good
1009 
1010  int32_t n_dim I # dimensions to use
1011  int32_t dat_ix I linear offset value
1012  int32_t * offset O offset array
1013 
1014 ********************************************************************/
1015 int32_t linear_to_offset( int32_t n_dim, int32_t dat_ix, int32_t *offset )
1016  {
1017  int32_t idim;
1018 
1019  for( idim = 0; idim < n_dim; idim++ )
1020  {
1021  /* the offset for that dim is the modulo 2 of the linear,
1022  then next time divide te linear offset by 2 to do the next dim offset */
1023  offset[idim] = dat_ix % 2;
1024  dat_ix /= 2;
1025  }
1026  return 0;
1027  }
int32_t init_dim(int32_t, int32_t, double *)
Definition: dim_mgr.cpp:172
an array had not been initialized Several spelling and grammar corrections were which is read from the appropriate MCF the above metadata values were hard coded A problem calculating the average background DN for SWIR bands when the moon is in the space view port was corrected The new algorithm used to calculate the average background DN for all reflective bands when the moon is in the space view port is now the same as the algorithm employed by the thermal bands For non SWIR changes in the averages are typically less than Also for non SWIR the black body DNs remain a backup in case the SV DNs are not available For SWIR the changes in computed averages were larger because the old which used the black body suffered from contamination by the micron leak As a consequence of the if SV DNs are not available for the SWIR the EV pixels will not be the granule time is used to identify the appropriate tables within the set given for one LUT the first two or last two tables respectively will be used for the interpolation If there is only one LUT in the set of it will be treated as a constant LUT The manner in which Earth View data is checked for saturation was changed Previously the raw Earth View DNs and Space View DNs were checked against the lookup table values contained in the table dn_sat The change made is to check the raw Earth and Space View DNs to be sure they are less than the maximum saturation value and to check the Space View subtracted Earth View dns against a set of values contained in the new lookup table dn_sat_ev The metadata configuration and ASSOCIATEDINSTRUMENTSHORTNAME from the MOD02HKM product The same metatdata with extensions and were removed from the MOD021KM and MOD02OBC products ASSOCIATEDSENSORSHORTNAME was set to MODIS in all products These changes are reflected in new File Specification which users may consult for exact the pow functions were eliminated in Emissive_Cal and Emissive bands replaced by more efficient code Other calculations throughout the code were also made more efficient Aside from a few round off there was no difference to the product The CPU time decreased by about for a day case and for a night case A minor bug in calculating the uncertainty index for emissive bands was corrected The frame index(0-based) was previously being used the frame number(1-based) should have been used. There were only a few minor changes to the uncertainty index(maximum of 1 digit). 3. Some inefficient arrays(Sigma_RVS_norm_sq) were eliminated and some code lines in Preprocess_L1A_Data were moved into Process_OBCEng_Emiss. There were no changes to the product. Required RAM was reduced by 20 MB. Now
void ** dat_ptrs
Definition: dim_mgr.hpp:70
int j
Definition: decode_rs.h:73
int status
Definition: l1_czcs_hdf.c:32
~dim_mgr()
Definition: dim_mgr.cpp:83
int32_t n_intervals
Definition: dim_mgr.hpp:24
These are used to scale the SD before writing it to the HDF4 file The default is and which means the product is not scaled at all Since the product is usually stored as a float inside of this is a way to write the float out as a integer l2prod min
double min
Definition: dim_mgr.hpp:14
double * wt
Definition: dim_mgr.hpp:68
#define NULL
Definition: decode_rs.h:63
valarray< int > ix_arr
Definition: dim_mgr.hpp:45
data_info_struc ** dat_info
Definition: dim_mgr.hpp:39
int32_t access_id
Definition: dim_mgr.hpp:33
int32_t prune(int32_t)
Definition: dim_mgr.cpp:821
valarray< int > ix_arr
Definition: dim_mgr.hpp:53
double * wt_pt
Definition: dim_mgr.hpp:69
int32_t * pt_status
Definition: dim_mgr.hpp:63
data_info_struc * dat_info
Definition: dim_mgr.hpp:54
void dump_last_interval()
Definition: dim_mgr.cpp:962
int32_t get_tot_blobs()
Definition: dim_mgr.cpp:123
int32_t nvals
Definition: dim_mgr.hpp:13
double * dim_coords
Definition: dim_mgr.hpp:18
int32_t get_hdim()
Definition: dim_mgr.cpp:999
void add_pts(int32_t)
Definition: dim_mgr.cpp:113
int32_t type
Definition: dim_mgr.hpp:17
int32_t * ix_arr
Definition: dim_mgr.hpp:27
pt_info_struc * mng_pt(double *, int32_t, int32_t *)
Definition: dim_mgr.cpp:246
Definition: dim_mgr.hpp:43
int32_t get_ndim()
Definition: dim_mgr.cpp:998
void update_new_data()
Definition: dim_mgr.cpp:202
int32_t interval_needs_data
Definition: dim_mgr.hpp:62
int32_t dat_status
Definition: dim_mgr.hpp:25
bool match(char *first, char *second)
Definition: l2qc_viirs.cpp:45
double max
Definition: dim_mgr.hpp:15
dim_mgr(int32_t)
Definition: dim_mgr.cpp:26
int32_t dat_filled
Definition: dim_mgr.hpp:37
void insert(std::unordered_set< std::type_index > &found_types)
Definition: Product.hpp:75
int32_t linear_to_offset(int32_t n_dim, int32_t dat_ix, int32_t *offset)
Definition: dim_mgr.cpp:1015
void dump_mgr(double *)
Definition: dim_mgr.cpp:897
int32_t nvals
data_t s[NROOTS]
Definition: decode_rs.h:75
l2prod offset
int32_t * pt_base_loc
Definition: dim_mgr.hpp:66
interval_struc * int_str
Definition: dim_mgr.hpp:46
int32_t purge()
Definition: dim_mgr.cpp:752
double delta
Definition: dim_mgr.hpp:16
int i
Definition: decode_rs.h:71
l2prod max