← Back to branch summary

~ubuntu-branches/ubuntu/saucy/darktable/saucy

« back to all changes in this revision

Viewing changes to src/iop/Permutohedral.h

  • Committer: Bazaar Package Importer
  • Author(s): David Bremner
  • Date: 2011-07-12 09:36:46 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20110712093646-yp9dbxan44dmw15h
Tags: 0.9-1
* New upstream release.
* Remove all patches now upstream; only patch for
  -Wno-error=unused-but-set-variable remains.
* Bump Standards-Version to 3.9.2 (no changes)

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/iop/Permutohedral.h
38
38
 * dimensional space.                                              *
39
39
 *                                                                 *
40
40
 *******************************************************************/
 
41
template <int KD, int VD>
41
42
class HashTablePermutohedral
42
43
{
43
44
public:
45
46
   *  kd_: the dimensionality of the position vectors on the hyperplane.
46
47
   *  vd_: the dimensionality of the value vectors
47
48
   */
48
 
  HashTablePermutohedral(int kd_, int vd_) : kd(kd_), vd(vd_)
 
49
  HashTablePermutohedral()
49
50
  {
50
51
    capacity = 1 << 15;
 
52
    capacity_bits = 0x7fff;
51
53
    filled = 0;
52
54
    entries = new Entry[capacity];
53
 
    keys = new short[kd*capacity/2];
54
 
    values = new float[vd*capacity/2];
55
 
    memset(values, 0, sizeof(float)*vd*capacity/2);
 
55
    keys = new short[KD*capacity/2];
 
56
    values = new float[VD*capacity/2];
 
57
    memset(values, 0, sizeof(float)*VD*capacity/2);
56
58
  }
57
59
 
58
60
  ~HashTablePermutohedral()
69
71
  }
70
72
 
71
73
  // Returns a pointer to the keys array.
72
 
  short *getKeys()
 
74
  const short *getKeys()
73
75
  {
74
76
    return keys;
75
77
  }
86
88
   *  create: a flag specifying whether an entry should be created,
87
89
   *          should an entry with the given key not found.
88
90
   */
89
 
  int lookupOffset(short *key, size_t h, bool create = true)
 
91
  int lookupOffset(const short *key, size_t h, bool create = true)
90
92
  {
91
93
 
92
94
    // Double hash table size if necessary
104
106
      {
105
107
        if (!create) return -1; // Return not found.
106
108
        // need to create an entry. Store the given key.
107
 
        for (int i = 0; i < kd; i++)
108
 
          keys[filled*kd+i] = key[i];
109
 
        e.keyIdx = filled*kd;
110
 
        e.valueIdx = filled*vd;
 
109
        for (int i = 0; i < KD; i++)
 
110
          keys[filled*KD+i] = key[i];
 
111
        e.keyIdx = filled*KD;
 
112
        e.valueIdx = filled*VD;
111
113
        entries[h] = e;
112
114
        filled++;
113
115
        return e.valueIdx;
115
117
 
116
118
      // check if the cell has a matching key
117
119
      bool match = true;
118
 
      for (int i = 0; i < kd && match; i++)
 
120
      for (int i = 0; i < KD && match; i++)
119
121
        match = keys[e.keyIdx+i] == key[i];
120
122
      if (match)
121
123
        return e.valueIdx;
130
132
   *        k : pointer to the key vector to be looked up.
131
133
   *   create : true if a non-existing key should be created.
132
134
   */
133
 
  float *lookup(short *k, bool create = true)
 
135
  float *lookup(const short *k, bool create = true)
134
136
  {
135
 
    size_t h = hash(k) % capacity;
 
137
    size_t h = hash(k) & capacity_bits;
136
138
    int offset = lookupOffset(k, h, create);
137
139
    if (offset < 0) return NULL;
138
140
    else return values + offset;
142
144
  size_t hash(const short *key)
143
145
  {
144
146
    size_t k = 0;
145
 
    for (int i = 0; i < kd; i++)
 
147
    for (int i = 0; i < KD; i++)
146
148
    {
147
149
      k += key[i];
148
150
      k *= 2531011;
158
160
 
159
161
    size_t oldCapacity = capacity;
160
162
    capacity *= 2;
 
163
    capacity_bits = (capacity_bits << 1) | 1;
161
164
 
162
165
    // Migrate the value vectors.
163
 
    float *newValues = new float[vd*capacity/2];
164
 
    memset(newValues, 0, sizeof(float)*vd*capacity/2);
165
 
    memcpy(newValues, values, sizeof(float)*vd*filled);
 
166
    float *newValues = new float[VD*capacity/2];
 
167
    memset(newValues, 0, sizeof(float)*VD*capacity/2);
 
168
    memcpy(newValues, values, sizeof(float)*VD*filled);
166
169
    delete[] values;
167
170
    values = newValues;
168
171
 
169
172
    // Migrate the key vectors.
170
 
    short *newKeys = new short[kd*capacity/2];
171
 
    memcpy(newKeys, keys, sizeof(short)*kd*filled);
 
173
    short *newKeys = new short[KD*capacity/2];
 
174
    memcpy(newKeys, keys, sizeof(short)*KD*filled);
172
175
    delete[] keys;
173
176
    keys = newKeys;
174
177
 
178
181
    for (size_t i = 0; i < oldCapacity; i++)
179
182
    {
180
183
      if (entries[i].keyIdx == -1) continue;
181
 
      size_t h = hash(keys + entries[i].keyIdx) % capacity;
 
184
      size_t h = hash(keys + entries[i].keyIdx) & capacity_bits;
182
185
      while (newEntries[h].keyIdx != -1)
183
186
      {
184
187
        h++;
202
205
  float *values;
203
206
  Entry *entries;
204
207
  size_t capacity, filled;
205
 
  int kd, vd;
 
208
  unsigned long capacity_bits;
206
209
};
207
210
 
208
211
/******************************************************************
211
214
 * PermutohedralLattice::filter(...) does all the work.           *
212
215
 *                                                                *
213
216
 ******************************************************************/
 
217
template <int D, int VD>
214
218
class PermutohedralLattice
215
219
{
216
220
public:
217
 
 
218
 
#if 0
219
 
  /* Filters given image against a reference image.
220
 
   *   im : image to be filtered.
221
 
   *  ref : reference image whose edges are to be respected.
222
 
   */
223
 
  static Image filter(Image im, Image ref)
224
 
  {
225
 
 
226
 
 
227
 
    // Create lattice
228
 
    PermutohedralLattice lattice(ref.channels, im.channels+1, im.width*im.height*im.frames);
229
 
 
230
 
    // Splat into the lattice
231
 
    //printf("Splatting...\n");
232
 
 
233
 
    float *imPtr = im(0, 0, 0);
234
 
    float *refPtr = ref(0, 0, 0);
235
 
    for (int t = 0; t < im.frames; t++)
236
 
    {
237
 
      for (int y = 0; y < im.height; y++)
238
 
      {
239
 
        for (int x = 0; x < im.width; x++)
240
 
        {
241
 
          lattice.splat(refPtr, imPtr);
242
 
          refPtr += ref.channels;
243
 
          imPtr += im.channels;
244
 
        }
245
 
      }
246
 
    }
247
 
 
248
 
    // Blur the lattice
249
 
    //printf("Blurring...\n");
250
 
    lattice.blur();
251
 
 
252
 
    // Slice from the lattice
253
 
    //printf("Slicing...\n");
254
 
 
255
 
    Image out(im.width, im.height, im.frames, im.channels);
256
 
 
257
 
    lattice.beginSlice();
258
 
    float *outPtr = out(0, 0, 0);
259
 
    for (int t = 0; t < im.frames; t++)
260
 
    {
261
 
      for (int y = 0; y < im.height; y++)
262
 
      {
263
 
        for (int x = 0; x < im.width; x++)
264
 
        {
265
 
          lattice.slice(outPtr);
266
 
          outPtr += out.channels;
267
 
        }
268
 
      }
269
 
    }
270
 
 
271
 
    return out;
272
 
  }
273
 
#endif
274
 
 
275
221
  /* Constructor
276
222
   *     d_ : dimensionality of key vectors
277
223
   *    vd_ : dimensionality of value vectors
278
224
   * nData_ : number of points in the input
279
225
   */
280
 
  PermutohedralLattice(int d_, int vd_, int nData_) :
281
 
    d(d_), vd(vd_), nData(nData_), hashTable(d_, vd_)
 
226
  PermutohedralLattice(int nData_, int nThreads_=1) :
 
227
    nData(nData_), nThreads(nThreads_)
282
228
  {
283
229
 
284
230
    // Allocate storage for various arrays
285
 
    elevated = new float[d+1];
286
 
    scaleFactor = new float[d];
 
231
    float *scaleFactorTmp = new float[D];
 
232
    int *canonicalTmp = new int[(D+1)*(D+1)];
287
233
 
288
 
    greedy = new short[d+1];
289
 
    rank = new char[d+1];
290
 
    barycentric = new float[d+2];
291
 
    replay = new ReplayEntry[nData*(d+1)];
292
 
    nReplay = 0;
293
 
    canonical = new short[(d+1)*(d+1)];
294
 
    key = new short[d+1];
 
234
    replay = new ReplayEntry[nData*(D+1)];
295
235
 
296
236
    // compute the coordinates of the canonical simplex, in which
297
237
    // the difference between a contained point and the zero
298
238
    // remainder vertex is always in ascending order. (See pg.4 of paper.)
299
 
    for (int i = 0; i <= d; i++)
 
239
    for (int i = 0; i <= D; i++)
300
240
    {
301
 
      for (int j = 0; j <= d-i; j++)
302
 
        canonical[i*(d+1)+j] = i;
303
 
      for (int j = d-i+1; j <= d; j++)
304
 
        canonical[i*(d+1)+j] = i - (d+1);
 
241
      for (int j = 0; j <= D-i; j++)
 
242
        canonicalTmp[i*(D+1)+j] = i;
 
243
      for (int j = D-i+1; j <= D; j++)
 
244
        canonicalTmp[i*(D+1)+j] = i - (D+1);
305
245
    }
 
246
    canonical = canonicalTmp;
306
247
 
307
248
    // Compute parts of the rotation matrix E. (See pg.4-5 of paper.)
308
 
    for (int i = 0; i < d; i++)
 
249
    for (int i = 0; i < D; i++)
309
250
    {
310
251
      // the diagonal entries for normalization
311
 
      scaleFactor[i] = 1.0f/(sqrtf((float)(i+1)*(i+2)));
 
252
      scaleFactorTmp[i] = 1.0f/(sqrtf((float)(i+1)*(i+2)));
312
253
 
313
254
      /* We presume that the user would like to do a Gaussian blur of standard deviation
314
255
       * 1 in each dimension (or a total variance of d, summed over dimensions.)
322
263
       *
323
264
       * So we need to scale the space by (d+1)sqrt(2/3).
324
265
       */
325
 
      scaleFactor[i] *= (d+1)*sqrtf(2.0/3);
 
266
      scaleFactorTmp[i] *= (D+1)*sqrtf(2.0/3);
326
267
    }
 
268
    scaleFactor = scaleFactorTmp;
 
269
 
 
270
    hashTables = new HashTablePermutohedral<D,VD>[nThreads];
327
271
  }
328
272
 
329
273
 
330
274
  ~PermutohedralLattice()
331
275
  {
332
276
    delete[] scaleFactor;
333
 
    delete[] elevated;
334
 
    delete[] greedy;
335
 
    delete[] rank;
336
 
    delete[] barycentric;
337
277
    delete[] replay;
338
278
    delete[] canonical;
339
 
    delete[] key;
 
279
    delete[] hashTables;
340
280
  }
341
281
 
342
282
 
343
283
  /* Performs splatting with given position and value vectors */
344
 
  void splat(float *position, float *value)
 
284
  void splat(float *position, float *value, int replay_index, int thread_index=0)
345
285
  {
 
286
    float elevated[D+1];
 
287
    int greedy[D+1];
 
288
    int rank[D+1];
 
289
    float barycentric[D+2];
 
290
    short key[D];
346
291
 
347
292
    // first rotate position into the (d+1)-dimensional hyperplane
348
 
    elevated[d] = -d*position[d-1]*scaleFactor[d-1];
349
 
    for (int i = d-1; i > 0; i--)
 
293
    elevated[D] = -D*position[D-1]*scaleFactor[D-1];
 
294
    for (int i = D-1; i > 0; i--)
350
295
      elevated[i] = (elevated[i+1] -
351
296
                     i*position[i-1]*scaleFactor[i-1] +
352
297
                     (i+2)*position[i]*scaleFactor[i]);
353
298
    elevated[0] = elevated[1] + 2*position[0]*scaleFactor[0];
354
299
 
355
300
    // prepare to find the closest lattice points
356
 
    float scale = 1.0f/(d+1);
357
 
    char * myrank = rank;
358
 
    short * mygreedy = greedy;
 
301
    float scale = 1.0f/(D+1);
359
302
 
360
303
    // greedily search for the closest zero-colored lattice point
361
304
    int sum = 0;
362
 
    for (int i = 0; i <= d; i++)
 
305
    for (int i = 0; i <= D; i++)
363
306
    {
364
307
      float v = elevated[i]*scale;
365
 
      float up = ceilf(v)*(d+1);
366
 
      float down = floorf(v)*(d+1);
367
 
 
368
 
      if (up - elevated[i] < elevated[i] - down) mygreedy[i] = (short)up;
369
 
      else mygreedy[i] = (short)down;
370
 
 
371
 
      sum += mygreedy[i];
 
308
      float up = ceilf(v)*(D+1);
 
309
      float down = floorf(v)*(D+1);
 
310
 
 
311
      if (up - elevated[i] < elevated[i] - down) greedy[i] = up;
 
312
      else greedy[i] = down;
 
313
 
 
314
      sum += greedy[i];
372
315
    }
373
 
    sum /= d+1;
 
316
    sum /= D+1;
374
317
 
375
318
    // rank differential to find the permutation between this simplex and the canonical one.
376
319
    // (See pg. 3-4 in paper.)
377
 
    memset(myrank, 0, sizeof(char)*(d+1));
378
 
    for (int i = 0; i < d; i++)
379
 
      for (int j = i+1; j <= d; j++)
380
 
        if (elevated[i] - mygreedy[i] < elevated[j] - mygreedy[j]) myrank[i]++;
381
 
        else myrank[j]++;
 
320
    memset(rank, 0, sizeof rank);
 
321
    for (int i = 0; i < D; i++)
 
322
      for (int j = i+1; j <= D; j++)
 
323
        if (elevated[i] - greedy[i] < elevated[j] - greedy[j]) rank[i]++;
 
324
        else rank[j]++;
382
325
 
383
326
    if (sum > 0)
384
327
    {
385
328
      // sum too large - the point is off the hyperplane.
386
329
      // need to bring down the ones with the smallest differential
387
 
      for (int i = 0; i <= d; i++)
 
330
      for (int i = 0; i <= D; i++)
388
331
      {
389
 
        if (myrank[i] >= d + 1 - sum)
 
332
        if (rank[i] >= D + 1 - sum)
390
333
        {
391
 
          mygreedy[i] -= d+1;
392
 
          myrank[i] += sum - (d+1);
 
334
          greedy[i] -= D+1;
 
335
          rank[i] += sum - (D+1);
393
336
        }
394
337
        else
395
 
          myrank[i] += sum;
 
338
          rank[i] += sum;
396
339
      }
397
340
    }
398
341
    else if (sum < 0)
399
342
    {
400
343
      // sum too small - the point is off the hyperplane
401
344
      // need to bring up the ones with largest differential
402
 
      for (int i = 0; i <= d; i++)
 
345
      for (int i = 0; i <= D; i++)
403
346
      {
404
 
        if (myrank[i] < -sum)
 
347
        if (rank[i] < -sum)
405
348
        {
406
 
          mygreedy[i] += d+1;
407
 
          myrank[i] += (d+1) + sum;
 
349
          greedy[i] += D+1;
 
350
          rank[i] += (D+1) + sum;
408
351
        }
409
352
        else
410
 
          myrank[i] += sum;
 
353
          rank[i] += sum;
411
354
      }
412
355
    }
413
356
 
414
357
    // Compute barycentric coordinates (See pg.10 of paper.)
415
 
    memset(barycentric, 0, sizeof(float)*(d+2));
416
 
    for (int i = 0; i <= d; i++)
 
358
    memset(barycentric, 0, sizeof barycentric);
 
359
    for (int i = 0; i <= D; i++)
417
360
    {
418
 
      barycentric[d-myrank[i]] += (elevated[i] - mygreedy[i]) * scale;
419
 
      barycentric[d+1-myrank[i]] -= (elevated[i] - mygreedy[i]) * scale;
 
361
      barycentric[D-rank[i]] += (elevated[i] - greedy[i]) * scale;
 
362
      barycentric[D+1-rank[i]] -= (elevated[i] - greedy[i]) * scale;
420
363
    }
421
 
    barycentric[0] += 1.0f + barycentric[d+1];
 
364
    barycentric[0] += 1.0f + barycentric[D+1];
422
365
 
423
366
    // Splat the value into each vertex of the simplex, with barycentric weights.
424
 
    for (int remainder = 0; remainder <= d; remainder++)
 
367
    for (int remainder = 0; remainder <= D; remainder++)
425
368
    {
426
369
      // Compute the location of the lattice point explicitly (all but the last coordinate - it's redundant because they sum to zero)
427
 
      for (int i = 0; i < d; i++)
428
 
        key[i] = mygreedy[i] + canonical[remainder*(d+1) + myrank[i]];
 
370
      for (int i = 0; i < D; i++)
 
371
        key[i] = greedy[i] + canonical[remainder*(D+1) + rank[i]];
429
372
 
430
373
      // Retrieve pointer to the value at this vertex.
431
 
      float * val = hashTable.lookup(key, true);
 
374
      float * val = hashTables[thread_index].lookup(key, true);
432
375
 
433
376
      // Accumulate values with barycentric weight.
434
 
      for (int i = 0; i < vd; i++)
 
377
      for (int i = 0; i < VD; i++)
435
378
        val[i] += barycentric[remainder]*value[i];
436
379
 
437
380
      // Record this interaction to use later when slicing
438
 
      replay[nReplay].offset = val - hashTable.getValues();
439
 
      replay[nReplay].weight = barycentric[remainder];
440
 
      nReplay++;
441
 
 
 
381
      replay[replay_index*(D+1)+remainder].table = thread_index;
 
382
      replay[replay_index*(D+1)+remainder].offset = val - hashTables[thread_index].getValues();
 
383
      replay[replay_index*(D+1)+remainder].weight = barycentric[remainder];
442
384
    }
443
385
  }
444
386
 
445
 
  // Prepare for slicing
446
 
  void beginSlice()
 
387
  /* Merge the multiple threads' hash tables into the totals. */
 
388
  void merge_splat_threads(void)
447
389
  {
448
 
    nReplay = 0;
 
390
    if (nThreads <= 1)
 
391
      return;
 
392
 
 
393
    /* Merge the multiple hash tables into one, creating an offset remap table. */
 
394
    int *offset_remap[nThreads];
 
395
    for (int i = 1; i < nThreads; i++)
 
396
    {
 
397
      const short *oldKeys = hashTables[i].getKeys();
 
398
      const float *oldVals = hashTables[i].getValues();
 
399
      const int filled = hashTables[i].size();
 
400
      offset_remap[i] = new int[filled];
 
401
      for (int j = 0; j < filled; j++)
 
402
      {
 
403
        float *val = hashTables[0].lookup(oldKeys+j*D, true);
 
404
        const float *oldVal = oldVals + j*VD;
 
405
        for (int k = 0; k < VD; k++)
 
406
          val[k] += oldVal[k];
 
407
        offset_remap[i][j] = val - hashTables[0].getValues();
 
408
      }
 
409
    }
 
410
 
 
411
    /* Rewrite the offsets in the replay structure from the above generated table. */
 
412
    for (int i = 0; i < nData*(D+1); i++)
 
413
      if (replay[i].table > 0)
 
414
        replay[i].offset = offset_remap[replay[i].table][replay[i].offset/VD];
 
415
 
 
416
    for (int i = 1; i < nThreads; i++)
 
417
      delete[] offset_remap[i];
449
418
  }
450
419
 
451
420
  /* Performs slicing out of position vectors. Note that the barycentric weights and the simplex
452
421
   * containing each position vector were calculated and stored in the splatting step.
453
422
   * We may reuse this to accelerate the algorithm. (See pg. 6 in paper.)
454
423
   */
455
 
  void slice(float *col)
 
424
  void slice(float *col, int replay_index)
456
425
  {
457
 
    float *base = hashTable.getValues();
458
 
    for (int j = 0; j < vd; j++) col[j] = 0;
459
 
    for (int i = 0; i <= d; i++)
 
426
    float *base = hashTables[0].getValues();
 
427
    for (int j = 0; j < VD; j++) col[j] = 0;
 
428
    for (int i = 0; i <= D; i++)
460
429
    {
461
 
      ReplayEntry r = replay[nReplay++];
462
 
      for (int j = 0; j < vd; j++)
 
430
      ReplayEntry r = replay[replay_index*(D+1)+i];
 
431
      for (int j = 0; j < VD; j++)
463
432
      {
464
433
        col[j] += r.weight*base[r.offset + j];
465
434
      }
470
439
  void blur()
471
440
  {
472
441
    // Prepare arrays
473
 
    short *neighbor1 = new short[d+1];
474
 
    short *neighbor2 = new short[d+1];
475
 
    float *newValue = new float[vd*hashTable.size()];
476
 
    float *oldValue = hashTable.getValues();
 
442
    float *newValue = new float[VD*hashTables[0].size()];
 
443
    float *oldValue = hashTables[0].getValues();
477
444
    float *hashTableBase = oldValue;
478
445
 
479
 
    float *zero = new float[vd];
480
 
    for (int k = 0; k < vd; k++) zero[k] = 0;
 
446
    float zero[VD];
 
447
    for (int k = 0; k < VD; k++) zero[k] = 0;
481
448
 
482
449
    // For each of d+1 axes,
483
 
    for (int j = 0; j <= d; j++)
 
450
    for (int j = 0; j <= D; j++)
484
451
    {
 
452
#ifdef _OPENMP
 
453
#pragma omp parallel for shared(j, oldValue, newValue, hashTableBase, zero)
 
454
#endif
485
455
      // For each vertex in the lattice,
486
 
      for (int i = 0; i < hashTable.size(); i++)   // blur point i in dimension j
 
456
      for (int i = 0; i < hashTables[0].size(); i++)   // blur point i in dimension j
487
457
      {
488
 
        short *key    = hashTable.getKeys() + i*(d); // keys to current vertex
489
 
        for (int k = 0; k < d; k++)
 
458
        const short *key    = hashTables[0].getKeys() + i*(D); // keys to current vertex
 
459
        short neighbor1[D+1];
 
460
        short neighbor2[D+1];
 
461
        for (int k = 0; k < D; k++)
490
462
        {
491
463
          neighbor1[k] = key[k] + 1;
492
464
          neighbor2[k] = key[k] - 1;
493
465
        }
494
 
        neighbor1[j] = key[j] - d;
495
 
        neighbor2[j] = key[j] + d; // keys to the neighbors along the given axis.
 
466
        neighbor1[j] = key[j] - D;
 
467
        neighbor2[j] = key[j] + D; // keys to the neighbors along the given axis.
496
468
 
497
 
        float *oldVal = oldValue + i*vd;
498
 
        float *newVal = newValue + i*vd;
 
469
        float *oldVal = oldValue + i*VD;
 
470
        float *newVal = newValue + i*VD;
499
471
 
500
472
        float *vm1, *vp1;
501
473
 
502
 
        vm1 = hashTable.lookup(neighbor1, false); // look up first neighbor
 
474
        vm1 = hashTables[0].lookup(neighbor1, false); // look up first neighbor
503
475
        if (vm1) vm1 = vm1 - hashTableBase + oldValue;
504
476
        else vm1 = zero;
505
477
 
506
 
        vp1 = hashTable.lookup(neighbor2, false); // look up second neighbor
 
478
        vp1 = hashTables[0].lookup(neighbor2, false); // look up second neighbor
507
479
        if (vp1) vp1 = vp1 - hashTableBase + oldValue;
508
480
        else vp1 = zero;
509
481
 
510
482
        // Mix values of the three vertices
511
 
        for (int k = 0; k < vd; k++)
 
483
        for (int k = 0; k < VD; k++)
512
484
          newVal[k] = (0.25f*vm1[k] + 0.5f*oldVal[k] + 0.25f*vp1[k]);
513
485
      }
514
486
      float *tmp = newValue;
520
492
    // depending where we ended up, we may have to copy data
521
493
    if (oldValue != hashTableBase)
522
494
    {
523
 
      memcpy(hashTableBase, oldValue, hashTable.size()*vd*sizeof(float));
 
495
      memcpy(hashTableBase, oldValue, hashTables[0].size()*VD*sizeof(float));
524
496
      delete[] oldValue;
525
497
    }
526
498
    else
527
499
    {
528
500
      delete[] newValue;
529
501
    }
530
 
 
531
 
    delete[] zero;
532
 
    delete[] neighbor1;
533
 
    delete[] neighbor2;
534
502
  }
535
503
 
536
504
private:
537
505
 
538
 
  int d, vd, nData;
539
 
  float *elevated, *scaleFactor, *barycentric;
540
 
  short *canonical;
541
 
  short *key;
 
506
  int nData;
 
507
  int nThreads;
 
508
  const float *scaleFactor;
 
509
  const int *canonical;
542
510
 
543
511
  // slicing is done by replaying splatting (ie storing the sparse matrix)
544
512
  struct ReplayEntry
545
513
  {
 
514
    int table;
546
515
    int offset;
547
516
    float weight;
548
517
  } *replay;
549
 
  int nReplay, nReplaySub;
550
518
 
551
 
public:
552
 
  char  *rank;
553
 
  short *greedy;
554
 
  HashTablePermutohedral hashTable;
 
519
  HashTablePermutohedral<D,VD> *hashTables;
555
520
};
556
521
 
557
522
#endif