1 //////////////////////////////////////////////////////////////////////
2 //
3 // Pixie
4 //
5 // Copyright � 1999 - 2003, Okan Arikan
6 //
7 // Contact: okan@cs.utexas.edu
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the GNU Lesser General Public
11 // License as published by the Free Software Foundation; either
12 // version 2.1 of the License, or (at your option) any later version.
13 //
14 // This library is distributed in the hope that it will be useful,
15 // but WITHOUT ANY WARRANTY; without even the implied warranty of
16 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 // Lesser General Public License for more details.
18 //
19 // You should have received a copy of the GNU Lesser General Public
20 // License along with this library; if not, write to the Free Software
21 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
22 //
23 ///////////////////////////////////////////////////////////////////////
24 ///////////////////////////////////////////////////////////////////////
25 //
26 // File : brickmap.cpp
27 // Classes : CBrickMap
28 // Description : Da implementation
29 //
30 ////////////////////////////////////////////////////////////////////////
31
32 #include "common/os.h"
33 #include "brickmap.h"
34 #include "random.h"
35 #include "error.h"
36 #include "pointCloud.h"
37 #include "memory.h"
38 #include "renderer.h"
39 #include "atomic.h"
40
41
42
43
44 // The static members of the CBrickMap class
45 CBrickMap *CBrickMap::brickMaps = NULL; // List of brickmaps in memory
46 int CBrickMap::referenceNumber = 0; // The last reference number
47 int CBrickMap::currentMemory = 0; // The currently used memory abount
48 int CBrickMap::maxMemory = 0; // The maximum memory for brickmaps
49 int CBrickMap::detailLevel = 2; // The detail level
50 int CBrickMap::drawType = 0; // Draw boxes
51 int CBrickMap::drawChannel = 0; // Which channel to draw
52
53
54 // convert brick to voxel
55 const float INV_BRICK_SIZE = 1.0f/ (float) BRICK_SIZE;
56 const float InvLog2 = 1.0f/log(2.0f);
57
58 ////////////////////////////////////////////////////////////////////////
59 // This number controls the density of the points in the finest level of the tree
60 // Larger it is, the smaller the leaf level voxels will be (compared to the point size)
61 //#define LEAF_FACTOR 0.4f
62
63 // the factor of 2 converts point radii to diameters, the factor of 1/8 from side to voxel size
64 const float LEAF_FACTOR = INV_BRICK_SIZE/2.0f;
65
66
67 ///////////////////////////////////////////////////////////////////////
68 // Function : intersect
69 // Description : Compute the volume of the intersection of the cube centered at P with side dP, with the cube centered at x,y,z with side d
70 // Return Value : the volume of the intersection
71 // Comments :
intersect(const float * P,float dP,float x,float y,float z,float d)72 static inline float intersect(const float *P,float dP,float x,float y,float z,float d) {
73 float tmin1,tmin2,tmax1,tmax2,tmin,tmax;
74 float w;
75
76 // X coordinate overlap
77 tmin1 = P[0] - dP;
78 tmax1 = P[0] + dP;
79 tmin2 = x - d;
80 tmax2 = x + d;
81 tmin = max(tmin1,tmin2);
82 tmax = min(tmax1,tmax2);
83 if (tmax <= tmin) return 0;
84 w = tmax - tmin;
85
86 // Y coordinate overlap
87 tmin1 = P[1] - dP;
88 tmax1 = P[1] + dP;
89 tmin2 = y - d;
90 tmax2 = y + d;
91 tmin = max(tmin1,tmin2);
92 tmax = min(tmax1,tmax2);
93 if (tmax <= tmin) return 0;
94 w *= tmax - tmin;
95
96 // Z coordinate overlap
97 tmin1 = P[2] - dP;
98 tmax1 = P[2] + dP;
99 tmin2 = z - d;
100 tmax2 = z + d;
101 tmin = max(tmin1,tmin2);
102 tmax = min(tmax1,tmax2);
103 if (tmax <= tmin) return 0;
104 w *= tmax - tmin;
105
106 return w;
107 }
108
109
110
111
112
113
114
115
116
117
118
119 ///////////////////////////////////////////////////////////////////////
120 // Class : CBrickMap
121 // Method : CBrickMap
122 // Description : Ctor
123 // Return Value : -
124 // Comments :
CBrickMap(FILE * in,const char * name,const float * from,const float * to)125 CBrickMap::CBrickMap(FILE *in,const char *name,const float *from,const float *to) : CTexture3d(name,from,to) {
126 int offset,i;
127
128 // Init the data
129 nextMap = brickMaps;
130 brickMaps = this;
131 normalThreshold = 0.7f;
132 file = in;
133 modifying = FALSE;
134 osCreateMutex(mutex);
135
136 // Read the header offset
137 fseek(file,-(long)sizeof(int),SEEK_END);
138 fread(&offset,1,sizeof(int),file);
139 fseek(file,offset,SEEK_SET);
140
141 // Read the class data
142 readChannels(file);
143
144 fread(bmin,1,sizeof(vector),file);
145 fread(bmax,1,sizeof(vector),file);
146 fread(center,1,sizeof(vector),file);
147 fread(&side,1,sizeof(float),file);
148 invSide = 1 / side;
149 fread(&maxDepth,1,sizeof(int),file);
150 fread(activeBricks,BRICK_HASHSIZE,sizeof(CBrickNode *),file);
151
152 // Read the permanent bricks nodes
153 for (i=0;i<BRICK_HASHSIZE;i++) {
154 if (activeBricks[i] != NULL) {
155 activeBricks[i] = NULL;
156
157 while(TRUE) {
158 CBrickNode *cNode = new CBrickNode;
159
160 fread(cNode,1,sizeof(CBrickNode),file);
161
162 assert(cNode->brick == NULL);
163 assert(cNode->fileIndex != -1);
164
165 if (cNode->next != NULL) {
166 cNode->next = activeBricks[i];
167 activeBricks[i] = cNode;
168 } else {
169 cNode->next = activeBricks[i];
170 activeBricks[i] = cNode;
171 break;
172 }
173 }
174 }
175 }
176 }
177
178
179 ///////////////////////////////////////////////////////////////////////
180 // Class : CBrickMap
181 // Method : CBrickMap
182 // Description : Ctor
183 // Return Value : -
184 // Comments : Use this contructor to compute from sctratch
CBrickMap(const char * name,const float * bmi,const float * bma,const float * from,const float * to,const float * toNDC,CChannel * ch,int nc,int md=10)185 CBrickMap::CBrickMap(const char *name,const float *bmi,const float *bma,const float *from,const float *to,const float *toNDC,CChannel *ch,int nc,int md = 10) : CTexture3d(name,from,to,toNDC,nc,ch) {
186 int i;
187
188 // Init the data
189 nextMap = brickMaps;
190 brickMaps = this;
191 normalThreshold = 0.7f;
192 file = NULL;
193 modifying = TRUE;
194 osCreateMutex(mutex);
195
196
197 // Compute the bounding cube
198 movvv(bmin,bmi);
199 movvv(bmax,bma);
200 subvv(bmax,bmin);
201 side = bmax[0];
202 side = max(side,bmax[1]);
203 side = max(side,bmax[2]);
204 invSide = 1 / side;
205 addvf(bmax,bmin,side);
206 addvv(center,bmin,bmax);
207 mulvf(center,0.5f);
208
209 maxDepth = md;
210 file = ropen(name,"wb+",fileBrickMap); // This is the file we will be writing to
211
212 // Initialize the hash table
213 for (i=0;i<BRICK_HASHSIZE;i++) activeBricks[i] = NULL;
214 }
215
216 ///////////////////////////////////////////////////////////////////////
217 // Class : CBrickMap
218 // Method : ~CBrickMap
219 // Description : Dtor
220 // Return Value : -
221 // Comments :
~CBrickMap()222 CBrickMap::~CBrickMap() {
223 int i;
224 CBrickNode *cNode;
225 CBrickMap *cMap,*pMap;
226
227 // Flush the memory
228 flushBrickMap(TRUE);
229
230 // Remove us from the list of bricks in memory
231 for (pMap=NULL,cMap=brickMaps;cMap!=NULL;pMap=cMap,cMap=cMap->nextMap) {
232 if (cMap == this) {
233 if (pMap == NULL) brickMaps = nextMap;
234 else pMap->nextMap = nextMap;
235 break;
236 }
237 }
238
239 // Free the hash table
240 for (i=0;i<BRICK_HASHSIZE;i++) {
241 while((cNode=activeBricks[i]) != NULL) {
242 activeBricks[i] = cNode->next;
243 assert(cNode->brick == NULL);
244 delete cNode;
245 }
246 }
247
248 // Close the file if not already have done so
249 if (file != NULL) fclose(file);
250
251 osDeleteMutex(mutex);
252 }
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271 //////////////////////////////////////////////////////////////////
272 //
273 // Here are some convenient macros for locating bricks/voxels
274
275
276
277
278
279
280
281
282 ///////////////////////////////////////////////////////////
283 // This macro iterates over the bricks that intersect the normalized point P
284 // ---> Preconditions:
285 // P = normalized point
286 // dP = lookup radius
287 // ---> Within the loop:
288 // x,y,z = the brick coordinates
289 #define forEachBrick(__depth) \
290 const float tmp = invSide*(float) (1 << __depth); \
291 int xs = (int) floor((P[0]-dP)*tmp); \
292 int ys = (int) floor((P[1]-dP)*tmp); \
293 int zs = (int) floor((P[2]-dP)*tmp); \
294 int xe = (int) floor((P[0]+dP)*tmp); \
295 int ye = (int) floor((P[1]+dP)*tmp); \
296 int ze = (int) floor((P[2]+dP)*tmp); \
297 int x,y,z; \
298 if (xs < 0) xs = 0; \
299 if (ys < 0) ys = 0; \
300 if (zs < 0) zs = 0; \
301 if (xe >= (1 << __depth)) xe = (1 << __depth) - 1; \
302 if (ye >= (1 << __depth)) ye = (1 << __depth) - 1; \
303 if (ze >= (1 << __depth)) ze = (1 << __depth) - 1; \
304 for (x=xs;x<=xe;x++) for (y=ys;y<=ye;y++) for (z=zs;z<=ze;z++) {
305
306 ///////////////////////////////////////////////////////////
307 // This macro iterates over the voxels that intersects the normalized point P
308 // ---> Preconditions:
309 // P = normalized point
310 // dP = lookup radius
311 // cBrick = the current brick
312 // ---> Within the loop:
313 // cWeight = the weight of the voxel
314 // cVoxel = the voxel
315 // cX,cY,cZ = the center of the voxel
316 #define forEachVoxel(__x,__y,__z,__depth) \
317 const float cSide = side / (float) (1 << __depth); \
318 const float xS = cSide*__x; \
319 const float yS = cSide*__y; \
320 const float zS = cSide*__z; \
321 const float dVoxel = cSide * INV_BRICK_SIZE; \
322 const float invDvoxel = 1 / dVoxel; \
323 char *cData = (char *) cBrick->voxels; \
324 int xvs = (int) floor(((P[0] - dP) - xS) * invDvoxel); \
325 int yvs = (int) floor(((P[1] - dP) - yS) * invDvoxel); \
326 int zvs = (int) floor(((P[2] - dP) - zS) * invDvoxel); \
327 int xve = (int) floor(((P[0] + dP) - xS) * invDvoxel); \
328 int yve = (int) floor(((P[1] + dP) - yS) * invDvoxel); \
329 int zve = (int) floor(((P[2] + dP) - zS) * invDvoxel); \
330 int xv,yv,zv; \
331 if (xvs < 0) xvs = 0; \
332 if (yvs < 0) yvs = 0; \
333 if (zvs < 0) zvs = 0; \
334 if (xve >= BRICK_SIZE) xve = BRICK_SIZE-1; \
335 if (yve >= BRICK_SIZE) yve = BRICK_SIZE-1; \
336 if (zve >= BRICK_SIZE) zve = BRICK_SIZE-1; \
337 for (xv=xvs;xv<=xve;xv++) for (yv=yvs;yv<=yve;yv++) for (zv=zvs;zv<=zve;zv++) { \
338 CVoxel *cVoxel = (CVoxel *) (cData + (zv*BRICK_SIZE*BRICK_SIZE + yv*BRICK_SIZE + xv)*(sizeof(CVoxel) + dataSize*sizeof(float))); \
339 const float cX = (xS + (xv + 0.5f)*dVoxel); \
340 const float cY = (yS + (yv + 0.5f)*dVoxel); \
341 const float cZ = (zS + (zv + 0.5f)*dVoxel); \
342 const float cWeight = intersect(P,dP,cX,cY,cZ,dVoxel*0.5f); \
343 if (cWeight == 0) continue;
344
345
346
347
348
349
350
351
352
353
354 ///////////////////////////////////////////////////////////////////////
355 // Class : CBrickMap
356 // Method : store
357 // Description : Add a point into the structure
358 // Return Value : -
359 // Comments :
store(const float * data,const float * cP,const float * cN,float dP)360 void CBrickMap::store(const float *data,const float *cP,const float *cN,float dP) {
361 dP *= dPscale;
362
363 int depth = (int) ceil(log(side*LEAF_FACTOR/dP)*InvLog2); // This is the depth we want to add
364 CBrick *cBrick;
365 CBrickNode *cNode;
366 vector P,N;
367
368 depth = min(max(depth,0),maxDepth);
369
370 // First, transform the point to world coordinate system
371 mulmp(P,to,cP);
372 assert(inBox(bmin,bmax,P));
373 subvv(P,bmin);
374 mulmn(N,from,cN);
375 if (dotvv(N,N) > 0) normalizev(N);
376
377 // Lock the structure
378 osLock(mutex);
379
380 // Iterate over the bricks we want
381 forEachBrick(depth)
382 int cDepth,cx,cy,cz;
383
384 // Forcefully create the bricks higher up in the hierarchy and add this point to the voxel data
385 for (cx=x,cy=y,cz=z,cDepth=depth;cDepth>=0;cx=cx>>1,cy=cy>>1,cz=cz>>1,cDepth--) {
386
387 // Create the brick
388 cBrick = findBrick(cx,cy,cz,cDepth,TRUE,&cNode);
389 assert(cBrick != NULL);
390
391 // Iterate over the voxels that intersect this one
392 forEachVoxel(cx,cy,cz,cDepth)
393
394 // Find the voxel we want to record
395 while(TRUE) {
396 const float tmp = dotvv(N,cVoxel->N);
397 // The small factor corrects for null normals / blank voxels not being matched
398 // numerical instability causes the first condition to fail when it should pass
399 if (((tmp*tmp + 1.0e-9f) >= (normalThreshold*normalThreshold*dotvv(cVoxel->N,cVoxel->N))) && (tmp*tmp >= 0)) {
400 break;
401 } else {
402 if (cVoxel->next == NULL) {
403 float *data;
404 int i;
405
406 cVoxel->next = (CVoxel *) new char[sizeof(CVoxel) + dataSize*sizeof(float)];
407 cVoxel = cVoxel->next;
408 data = (float *) (cVoxel+1);
409 initv(cVoxel->N,0);
410 cVoxel->weight = 0;
411 cVoxel->next = NULL;
412 for (i=0;i<dataSize;i++) data[i] = 0;
413
414 currentMemory += sizeof(CVoxel) + dataSize*sizeof(float);
415
416 // Mark the brick as needing new storage
417 // Note: we will need to compact the map afterwards
418 cNode->fileIndex = -1;
419 break;
420 } else {
421 cVoxel = cVoxel->next;
422 }
423 }
424 }
425
426 float *dest = (float *) (cVoxel+1);
427 cVoxel->N[0] += N[0]*cWeight;
428 cVoxel->N[1] += N[1]*cWeight;
429 cVoxel->N[2] += N[2]*cWeight;
430 for (int j=0;j<dataSize;j++) dest[j] += data[j]*cWeight;
431 cVoxel->weight += cWeight;
432 }
433 }
434 }
435
436 // Release the structure
437 osUnlock(mutex);
438 }
439
440
441
442 ///////////////////////////////////////////////////////////////////////
443 // Class : CBrickMap
444 // Method : lookup
445 // Description : Lookup data
446 // Return Value : -
447 // Comments :
448 void CBrickMap::lookup(float *data,const float *cP,const float *cN,float dP) {
449 dP *= dPscale;
450
451 float depthf = log(side*LEAF_FACTOR/dP)*InvLog2;
452 int depth = (int) floor(depthf);
453 float t = depthf - depth;
454 float *data0 = (float *) alloca(dataSize*2*sizeof(float));
455 float *data1 = data0 + dataSize;
456 float normalFactor= 1.0f;
457 vector P,N;
458 int i;
459
460 // First, transform the point to world coordinate system
461 mulmp(P,to,cP);
462 //assert(inBox(bmin,bmax,P));
463 subvv(P,bmin);
464 mulmn(N,from,cN);
465
466 if (dotvv(N,N) > 0) normalizev(N);
467 else normalFactor = 0.0f;
468
469 if (depth < 0) {
470 depth = 0;
471 depthf = 0;
472 }
473
474 // Perform the lookup
475 osLock(mutex);
476 atomicIncrement(&stats.numBrickmapLookups);
477 atomicIncrement(&stats.numBrickmapLookups);
478 lookup(P,N,dP,data0,depth,normalFactor);
479 lookup(P,N,dP,data1,depth+1,normalFactor);
480 osUnlock(mutex);
481
482 for (i=0;i<dataSize;i++) data[i] = data0[i]*(1-t) + data1[i]*t;
483 }
484
485
486
487 ///////////////////////////////////////////////////////////////////////
488 // Class : CBrickMap
489 // Method : lookup
490 // Description : Lookup a particular depth
491 // Return Value : -
492 // Comments :
493 void CBrickMap::lookup(const float *P,const float *N,float dP,float *data,int depth,float normalFactor) {
494 CBrick *cBrick;
495 float totalWeight = 0;
496 int i;
497
498 // Clear the data
499 for (i=0;i<dataSize;i++) data[i] = 0;
500
501 // Find the brick we want to look at
502 forEachBrick(depth)
503 int cDepth,cx,cy,cz;
504
505 // iterate all levels until we hit a valid sample
506 for (cx=x,cy=y,cz=z,cDepth=depth;cDepth>=0;cx=cx>>1,cy=cy>>1,cz=cz>>1,cDepth--) {
507
508 // Get the current brick
509 if ((cBrick = findBrick(cx,cy,cz,cDepth,FALSE,NULL)) != NULL) {
510 forEachVoxel(cx,cy,cz,cDepth)
511
512 // Find the voxel with the closest normal
513 for (;cVoxel!=NULL;cVoxel=cVoxel->next) {
514 const float weight = cWeight*cVoxel->weight*(normalFactor*dotvv(cVoxel->N,N) + (1.0f-normalFactor));
515
516 if (weight > 0) {
517 int j;
518 const float *src = (float *) (cVoxel+1);
519
520 for (j=0;j<dataSize;j++) data[j] += src[j]*weight;
521 totalWeight += weight;
522 }
523 }
524 }
525 }
526
527 // If we hit anything, we're done
528 if(totalWeight > 0) break;
529 }
530 }
531
532 // Normalize the data
533 if (totalWeight > 0) {
534 totalWeight = 1/totalWeight;
535 for (i=0;i<dataSize;i++) data[i] *= totalWeight;
536 }
537 }
538
539 ///////////////////////////////////////////////////////////////////////
540 // Class : CBrickMap
541 // Method : finalize
542 // Description : Finalize the creation of the brickmap
543 // Return Value : -
544 // Comments : this relies on all levels lower than a finer one being present
545 void CBrickMap::finalize() {
546 int *stack = (int *) alloca(maxDepth*8*5*sizeof(int));
547 int *stackBase = stack;
548 int headerOffset;
549 CBrickNode *cNode;
550 int i;
551
552 *stack++ = 0;
553 *stack++ = 0;
554 *stack++ = 0;
555 *stack++ = 0;
556 while(stack > stackBase) {
557 int x,y,z,depth;
558 CBrick *cBrick;
559
560 depth = *(--stack);
561 z = *(--stack);
562 y = *(--stack);
563 x = *(--stack);
564
565 if ((cBrick=findBrick(x,y,z,depth,FALSE,NULL)) != NULL) {
566 CVoxel *cVoxel;
567 int i;
568
569 // Make sure we iterate over the children
570 #define push(__x,__y,__z,__depth) *stack++ = __x; *stack++ = __y; *stack++ = __z; *stack++ = __depth;
571 push(2*x, 2*y, 2*z, depth+1);
572 push(2*x+1, 2*y, 2*z, depth+1);
573 push(2*x, 2*y+1, 2*z, depth+1);
574 push(2*x+1, 2*y+1, 2*z, depth+1);
575 push(2*x, 2*y, 2*z+1, depth+1);
576 push(2*x+1, 2*y, 2*z+1, depth+1);
577 push(2*x, 2*y+1, 2*z+1, depth+1);
578 push(2*x+1, 2*y+1, 2*z+1, depth+1);
579 #undef push
580
581 // Normalize the voxel data
582 for (cVoxel=cBrick->voxels,i=BRICK_SIZE*BRICK_SIZE*BRICK_SIZE;i>0;i--) {
583 float *vdata = (float*) (cVoxel+1);
584
585 // Deal with normalizing incoherent normals data
586 while(TRUE) {
587 float *data = (float *) (cVoxel+1);
588 if (cVoxel->weight > 0) {
589 const float invWeight = 1 / cVoxel->weight;
590 int j;
591
592 if (dotvv(cVoxel->N,cVoxel->N) > 0) normalizev(cVoxel->N);
593 for (j=0;j<dataSize;j++) data[j] *= invWeight;
594 cVoxel->weight = 1;
595 }
596
597 if (cVoxel->next != NULL) {
598 cVoxel = cVoxel->next;
599 } else {
600 break;
601 }
602 }
603
604 cVoxel = (CVoxel*) (vdata + dataSize);
605 }
606 }
607 }
608
609 // Flush all the bricks to disk
610 flushBrickMap(TRUE);
611
612 // Save the current file position. This is the file index
613 fseek(file,0,SEEK_END);
614 headerOffset = ftell(file);
615
616 // Write the class data here
617 writeChannels(file);
618
619 // Write the class data
620 fwrite(bmin,sizeof(vector),1,file);
621 fwrite(bmax,sizeof(vector),1,file);
622 fwrite(center,sizeof(vector),1,file);
623 fwrite(&side,sizeof(float),1,file);
624 fwrite(&maxDepth,sizeof(int),1,file);
625 fwrite(activeBricks,sizeof(CBrickNode *),BRICK_HASHSIZE,file);
626 for (i=0;i<BRICK_HASHSIZE;i++) {
627 for (cNode=activeBricks[i];cNode!=NULL;cNode=cNode->next) {
628
629 // Make sure the node is written to disk
630 assert(cNode->brick == NULL);
631 assert(cNode->fileIndex != -1);
632
633 fwrite(cNode,sizeof(CBrickNode),1,file);
634 }
635 }
636
637 // Write the file header at the beginning
638 fwrite(&headerOffset,sizeof(int),1,file);
639
640 // Mark the map as non-modifying, meaning
641 // we will no longer page out nodes
642 // this provides a big speed increase when
643 // compressing
644 modifying = FALSE;
645 }
646
647
648
649
650
651
652
653
654
655
656
657
658
659 ///////////////////////////////////////////////////////////////////////
660 // Class : CBrickMap
661 // Method : newBrick
662 // Description : Create a new brick
663 // Return Value : -
664 // Comments :
665 CBrickMap::CBrick *CBrickMap::newBrick(int clear) {
666 CBrick *cBrick;
667
668 // If we're using too much memory, swap some bricks out
669 if (currentMemory > maxMemory) flushBrickMap(FALSE);
670
671 // Allocate the brick
672 cBrick = (CBrick *) new char[sizeof(CBrick) + (sizeof(CVoxel) + dataSize*sizeof(float))*(BRICK_SIZE*BRICK_SIZE*BRICK_SIZE)];
673 cBrick->voxels = (CVoxel *) (cBrick + 1);
674
675 // Update the used memory
676 currentMemory += sizeof(CBrick) + (sizeof(CVoxel) + dataSize*sizeof(float))*(BRICK_SIZE*BRICK_SIZE*BRICK_SIZE);
677
678 if (clear) {
679 CVoxel *cVoxel;
680 int i,j;
681
682 // Clear the voxels
683 for (cVoxel=cBrick->voxels,i=0;i<(BRICK_SIZE*BRICK_SIZE*BRICK_SIZE);i++) {
684 float *data;
685 initv(cVoxel->N,0);
686 cVoxel->weight = 0;
687 cVoxel->next = NULL;
688 data = (float *) (cVoxel+1);
689 for (j=0;j<dataSize;j++) data[j] = 0;
690 cVoxel = (CVoxel *) ((char *) cVoxel + sizeof(CVoxel) + dataSize*sizeof(float));
691 }
692 }
693
694 return cBrick;
695 }
696
697 ///////////////////////////////////////////////////////////////////////
698 // Class : CBrickMap
699 // Method : loadBrick
700 // Description : Load a brick
701 // Return Value : -
702 // Comments :
703 CBrickMap::CBrick *CBrickMap::loadBrick(int fileIndex) {
704 CBrick *cBrick = newBrick(FALSE);
705 CVoxel *cVoxel,*tVoxel;
706 int i,j;
707
708 atomicIncrement(&stats.numBrickmapCachePageins);
709
710 // Seek to the right position in file
711 if (file == NULL) file = ropen(name,"w+",fileBrickMap);
712 fseek(file,fileIndex,SEEK_SET);
713
714 uint32_t bs[BRICK_PRESENCE_LONGS];
715 uint32_t b;
716
717 // work out which top-level voxels are present
718 fread(bs,sizeof(uint32_t)*BRICK_PRESENCE_LONGS,1,file);
719
720 // read those that are
721 for(i=0,cVoxel=cBrick->voxels;i<BRICK_PRESENCE_LONGS;i++){
722 b = bs[i];
723
724 for(j=BRICK_VOXEL_BATCH;j>0;j--) {
725 float *vdata = (float*) (cVoxel + 1);
726
727 if (b & 0x80000000L) {
728 fread(cVoxel,sizeof(CVoxel) + sizeof(float)*dataSize,1,file);
729
730 if (cVoxel->next != NULL) {
731 cVoxel->next = NULL;
732
733 while (TRUE) {
734 tVoxel = (CVoxel*) new char[sizeof(CVoxel) + dataSize*sizeof(float)];
735 currentMemory += sizeof(CVoxel) + dataSize*sizeof(float);
736
737 fread(tVoxel,sizeof(CVoxel) + sizeof(float)*dataSize,1,file);
738
739 if (tVoxel->next != NULL) {
740 tVoxel->next = cVoxel->next;
741 cVoxel->next = tVoxel;
742 } else {
743 tVoxel->next = cVoxel->next;
744 cVoxel->next = tVoxel;
745 break;
746 }
747 }
748 }
749 } else {
750 // initialize to null any that are not
751 cVoxel->weight = 0;
752 cVoxel->next = NULL;
753 initv(cVoxel->N,0);
754 }
755
756 b = b<<1;
757
758 cVoxel = (CVoxel*) (vdata + dataSize);
759 }
760 }
761
762 if (currentMemory > stats.brickmapPeakMem) stats.brickmapPeakMem = currentMemory;
763
764 // Return the brick
765 return cBrick;
766 }
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782 ///////////////////////////////////////////////////////////////////////
783 // Class : CBrickMap
784 // Method : brickMapCompact
785 // Description : This function is called to re-claim some of the memory from the brickmap
786 // Return Value : -
787 // Comments :
788 void CBrickMap::compact(const char *outFileName,float maxVariation) {
789 int numNodes;
790 CBrickNode *cNode;
791 CVoxel *cVoxel,*tVoxel;
792 CBrick *cBrick;
793 int i,j,k,vCnt,nullCnt,numCulled;
794
795 FILE *outfile = ropen(outFileName,"wb+",fileBrickMap);
796
797 // Use our own temp memory so we don't get mutex-reentrancy locking issues
798 CMemPage *tempMemory = NULL;
799 memoryInit(tempMemory);
800 memBegin(tempMemory);
801
802 CVoxel *tempVoxel = (CVoxel*) ralloc(sizeof(CVoxel) + dataSize*sizeof(float),tempMemory);
803 CBrickNode **newHash = (CBrickNode**) ralloc(BRICK_HASHSIZE*sizeof(CBrickNode*),tempMemory);
804 float *dataMean = (float*) ralloc(2*dataSize*sizeof(float),tempMemory);
805 float *dataVar = dataMean + dataSize;
806
807 // Initialize the hash
808 for (i=0;i<BRICK_HASHSIZE;i++) newHash[i] = NULL;
809
810 // Collect the loaded bricks into an array
811 numCulled = 0;
812 numNodes = 0;
813 nullCnt = 0;
814 for (i=0;i<BRICK_HASHSIZE;i++) {
815
816 for (cNode=activeBricks[i];cNode!=NULL;cNode=cNode->next) {
817
818 // Make sure we have the data
819 if (cNode->brick == NULL) {
820 // Get the thing resident
821 cNode->brick = loadBrick(cNode->fileIndex);
822 cNode->brick->referenceNumber = referenceNumber;
823 }
824 cBrick = cNode->brick;
825
826 // Calculate variance
827
828 for (j=0;j<dataSize;j++) dataMean[j] = dataVar[j] = 0;
829
830 vCnt = 0;
831
832 for (cVoxel=cBrick->voxels,k=BRICK_SIZE*BRICK_SIZE*BRICK_SIZE;k>0;k--) {
833 float *vdata = (float*) (cVoxel+1);
834
835 // Deal with normalizing incoherent normals data
836 while(TRUE) {
837 float *data = (float *) (cVoxel+1);
838
839 if(cVoxel->weight >0) {
840 for (j=0;j<dataSize;j++) dataMean[j] += data[j];
841 vCnt++;
842 } else nullCnt++;
843
844 if (cVoxel->next != NULL) {
845 cVoxel = cVoxel->next;
846 } else {
847 break;
848 }
849 }
850
851 cVoxel = (CVoxel*) (vdata + dataSize);
852 }
853
854 // Skip if we have no data in this brick
855 if (vCnt == 0) {
856 numCulled++;
857 continue;
858 }
859
860 float invCnt = 1.0f/(float)vCnt;
861 for (j=0;j<dataSize;j++) dataMean[j] *= invCnt;
862
863 for (cVoxel=cBrick->voxels,k=BRICK_SIZE*BRICK_SIZE*BRICK_SIZE;k>0;k--) {
864 float *vdata = (float*) (cVoxel+1);
865
866 // Deal with normalizing incoherent normals data
867 while(TRUE) {
868 float *data = (float *) (cVoxel+1);
869
870 if(cVoxel->weight >0) {
871 for (j=0;j<dataSize;j++) {
872 float d = (data[j]-dataMean[j]);
873 dataVar[j] += d*d;
874 }
875 }
876
877 if (cVoxel->next != NULL) {
878 cVoxel = cVoxel->next;
879 } else {
880 break;
881 }
882 }
883
884 cVoxel = (CVoxel*) (vdata + dataSize);
885 }
886
887 float maxVar = 0;
888 for (j=0;j<dataSize;j++) {
889 dataVar[j] *= invCnt;
890 dataVar[j] = sqrtf(dataVar[j]);
891 dataVar[j] /= dataMean[j];
892 if (dataVar[j] > maxVar) maxVar = dataVar[j];
893 }
894
895 // Do not write this brick if variation too low
896 if (maxVar < maxVariation && cNode->d > 0) {
897 numCulled++;
898 continue;
899 }
900
901
902 CBrickNode *tNode = (CBrickNode*) ralloc(sizeof(CBrickNode),tempMemory);
903
904 // Initialize temporary node data over
905 *tNode = *cNode;
906 tNode->next = newHash[i];
907 newHash[i] = tNode;
908 tNode->brick = NULL;
909 numNodes++;
910
911 // write it out to a new location
912 fseek(outfile,0,SEEK_END);
913 tNode->fileIndex = ftell(outfile);
914
915 uint32_t bs[BRICK_PRESENCE_LONGS];
916 uint32_t b;
917
918 // Work out for each voxel if there is anything at all
919 for (k=0,cVoxel=cBrick->voxels;k<BRICK_PRESENCE_LONGS;k++) {
920 b = 0;
921
922 for (j=BRICK_VOXEL_BATCH;j>0;j--) {
923 float *vdata = (float*) (cVoxel + 1);
924
925 b = b<<1;
926
927 while (cVoxel != NULL) {
928 if (cVoxel->weight > 0){
929 b |= 1;
930 }
931 cVoxel = cVoxel->next;
932 }
933
934 cVoxel = (CVoxel*) (vdata + dataSize);
935 }
936
937 bs[k] = b;
938 }
939
940 // Write the voxel-presence bits
941 fwrite(bs,sizeof(uint32_t)*BRICK_PRESENCE_LONGS,1,outfile);
942
943 // Write each voxel which exists
944 for(j=BRICK_SIZE*BRICK_SIZE*BRICK_SIZE,cVoxel=cBrick->voxels;j>0;j--) {
945 float *vdata = (float*) (cVoxel + 1);
946
947 int skippedLast = FALSE;
948 tVoxel = NULL;
949
950 while (cVoxel != NULL) {
951
952 skippedLast = FALSE;
953
954 if (cVoxel->weight > 0) {
955 fwrite(cVoxel,sizeof(CVoxel)+sizeof(float)*dataSize,1,outfile);
956 tVoxel = cVoxel;
957 } else {
958 skippedLast = TRUE;
959 }
960
961 cVoxel = cVoxel->next;
962 }
963
964 if ((skippedLast == TRUE) && (tVoxel != NULL)) {
965 // That last voxel we wrote had a bad next value...
966 fseek(outfile,-(long)(sizeof(CVoxel)+dataSize*sizeof(float)),SEEK_CUR);
967 memcpy(tempVoxel,tVoxel,sizeof(CVoxel)+sizeof(float)*dataSize);
968 tempVoxel->next = NULL;
969 // only resave the voxel header
970 fwrite(tempVoxel,sizeof(CVoxel),1,outfile);
971 // but keep the file write pointer correct
972 fseek(outfile,sizeof(float)*dataSize,SEEK_CUR);
973 }
974
975 cVoxel = (CVoxel*) (vdata + dataSize);
976 }
977 }
978 }
979 //fprintf(stderr,"%d bricks culled. %d null voxels not written\n",numCulled,nullCnt);
980
981 // Write out temporary node hash
982 // Save the current file position. This is the file index
983 fseek(outfile,0,SEEK_END);
984 int headerOffset = ftell(outfile);
985
986 // Write the class data here
987 writeChannels(outfile);
988
989 fwrite(bmin,sizeof(vector),1,outfile);
990 fwrite(bmax,sizeof(vector),1,outfile);
991 fwrite(center,sizeof(vector),1,outfile);
992 fwrite(&side,sizeof(float),1,outfile);
993 fwrite(&maxDepth,sizeof(int),1,outfile);
994 fwrite(newHash,sizeof(CBrickNode *),BRICK_HASHSIZE,outfile);
995 for (i=0;i<BRICK_HASHSIZE;i++) {
996 for (cNode=newHash[i];cNode!=NULL;cNode=cNode->next) {
997
998 // Make sure the node is written to disk
999 assert(cNode->brick == NULL);
1000 assert(cNode->fileIndex != -1);
1001
1002 fwrite(cNode,sizeof(CBrickNode),1,outfile);
1003 }
1004 }
1005
1006 // Write the position of the file header right at the end
1007 fwrite(&headerOffset,sizeof(int),1,outfile);
1008
1009 fclose(outfile);
1010
1011 memEnd(tempMemory);
1012 memoryTini(tempMemory);
1013 }
1014
1015
1016
1017
1018 ///////////////////////////////////////////////////////////////////////
1019 // Class : CBrickMap
1020 // Method : draw
1021 // Description : Draw the brickmap
1022 // Return Value : -
1023 // Comments :
1024 void CBrickMap::draw() {
1025 float P[chunkSize*3];
1026 float C[chunkSize*3];
1027 float N[chunkSize*3];
1028 float R[chunkSize];
1029 int j;
1030 int sampleStart = channels[drawChannel].sampleStart;
1031 int numSamples = channels[drawChannel].numSamples;
1032 float *cP = P;
1033 float *cC = C;
1034 float *cN = N;
1035 float *cR = R;
1036 int level = min(max(0,detailLevel),maxDepth);
1037 int nb = 1 << level;
1038 const float sqrt2 = sqrtf(0.5f);
1039 float cubePoints[] = { 0, 0, 0,
1040 1, 0, 0,
1041 1, 0, 1,
1042 0, 0, 1,
1043
1044 0, 1, 0,
1045 1, 1, 0,
1046 1, 1, 1,
1047 0, 1, 1,
1048
1049 0, 0, 0,
1050 1, 0, 0,
1051 1, 1, 0,
1052 0, 1, 0,
1053
1054 0, 0, 1,
1055 1, 0, 1,
1056 1, 1, 1,
1057 0, 1, 1,
1058
1059 0, 0, 0,
1060 0, 1, 0,
1061 0, 1, 1,
1062 0, 0, 1,
1063
1064 1, 0, 0,
1065 1, 1, 0,
1066 1, 1, 1,
1067 1, 0, 1
1068 };
1069
1070 if (0) {
1071 // draw bounding box lines
1072 j = chunkSize;
1073 float *pts = cubePoints;
1074 for(int k =0; k<6; k++) {
1075
1076 #define emitLn(i) \
1077 if (j == 0) { \
1078 drawLines(chunkSize,P,C); \
1079 cP = P; \
1080 cC = C; \
1081 j = chunkSize; \
1082 } \
1083 mulvf(cP,pts+ i*3,side); \
1084 initv(cC,0); \
1085 cP += 3; \
1086 cC += 3; \
1087 j--;
1088
1089 emitLn(0);
1090 emitLn(1);
1091 emitLn(1);
1092 emitLn(2);
1093 emitLn(2);
1094 emitLn(3);
1095 emitLn(3);
1096 emitLn(0);
1097
1098 pts += 12;
1099 }
1100 if (j != chunkSize) {
1101 drawLines(chunkSize-j,P,C);
1102 cP = P;
1103 cC = C;
1104 j = chunkSize;
1105 }
1106
1107 {
1108 const float hd = side/2.0f;
1109 const float d = side+side/2.0f;
1110 float tmp[9] = { d,hd,hd, hd,d,hd, hd,hd,d};
1111 float tmpc[9] = { 0 };
1112
1113 drawPoints(3,tmp,tmpc);
1114 }
1115 }
1116
1117 // For each brick at this level
1118 j = chunkSize;
1119 for (int xe=0;xe<nb;xe++) for (int ye=0;ye<nb;ye++) for (int ze=0;ze<nb;ze++) {
1120 float sz = side/(float) nb;
1121 int x=xe,y=ye,z=ze;
1122 CBrickMap::CBrick *bk = findBrick(x,y,z,level,false,NULL);
1123
1124 if (bk == NULL) continue;
1125
1126 CBrickMap::CVoxel *vx = bk->voxels;
1127
1128 float *DD = (float*)((char*)vx + sizeof(CBrickMap::CVoxel));
1129
1130 // For each voxel
1131 for(int zi=0;zi<BRICK_SIZE;zi++) for(int yi=0;yi<BRICK_SIZE;yi++) for(int xi=0;xi<BRICK_SIZE;xi++) {
1132 vector cent,Ctmp;
1133
1134 initv(cent,x*sz + xi*sz*INV_BRICK_SIZE,y*sz + yi*sz*INV_BRICK_SIZE,z*sz + zi*sz*INV_BRICK_SIZE);
1135
1136 // Save values before we update
1137 float *DDs = DD + sampleStart;
1138 if (numSamples == 1) {
1139 initv(Ctmp,DDs[0]);
1140 DDs = Ctmp;
1141 } else if (numSamples == 2) {
1142 initv(Ctmp,DDs[0],DDs[1],0);
1143 DDs = Ctmp;
1144 }
1145 float wt = vx->weight;
1146 float *norm = vx->N;
1147
1148 // Update for next iteration, incase we skip
1149 vx = (CBrickMap::CVoxel*)((char*)vx + sizeof(float)*dataSize + sizeof(CBrickMap::CVoxel));
1150 DD = (float*)((char*) DD + sizeof(float)*dataSize + sizeof(CBrickMap::CVoxel));
1151
1152 if (wt <= C_EPSILON) continue;
1153
1154 if (drawType == 0) {
1155
1156 float *pts = cubePoints;
1157 for(int k =0; k<6; k++) {
1158 vector tmp;
1159
1160 #define emitPt(i) \
1161 if (j == 0) { \
1162 drawTriangles(chunkSize,P,C); \
1163 cP = P; \
1164 cC = C; \
1165 j = chunkSize; \
1166 } \
1167 mulvf(tmp,pts+ i*3,sz/8.0f); \
1168 addvv(tmp,cent); \
1169 movvv(cP,tmp); \
1170 movvv(cC,DDs); \
1171 cP += 3; \
1172 cC += 3; \
1173 j--;
1174
1175 emitPt(0);
1176 emitPt(1);
1177 emitPt(2);
1178 emitPt(2);
1179 emitPt(3);
1180 emitPt(0);
1181
1182 pts += 12;
1183 }
1184 }
1185 else if (drawType == 1) {
1186 if (j == 0) {
1187 drawDisks(chunkSize,P,R,N,C);
1188 cP = P;
1189 cC = C;
1190 cN = N;
1191 cR = R;
1192 j = chunkSize;
1193 }
1194
1195 movvv(cP,cent);
1196 movvv(cC,DDs);
1197 movvv(cN,norm);
1198 cR[0] = sqrt2*sz*INV_BRICK_SIZE;
1199
1200 cP += 3;
1201 cC += 3;
1202 cN += 3;
1203 cR += 1;
1204 j--;
1205 } else if (drawType == 2) {
1206 if (j == 0) {
1207 drawPoints(chunkSize,P,C);
1208 cP = P;
1209 cC = C;
1210 j = chunkSize;
1211 }
1212
1213 movvv(cP,cent);
1214 movvv(cC,DDs);
1215
1216 cP += 3;
1217 cC += 3;
1218 j--;
1219 }
1220 }
1221 }
1222
1223 if (j != chunkSize) {
1224 if (drawType == 0) drawTriangles(chunkSize-j,P,C);
1225 else if (drawType == 1) drawDisks(chunkSize-j,P,R,N,C);
1226 else drawPoints(chunkSize-j,P,C);
1227 }
1228 }
1229
1230
1231 ///////////////////////////////////////////////////////////////////////
1232 // Class : CBrickMap
1233 // Method : keyDown
1234 // Description : handle interface keys
1235 // Return Value : -
1236 // Comments :
1237 int CBrickMap::keyDown(int key) {
1238 if ((key == 'M') || (key == 'm')) {
1239 detailLevel++;
1240 printf("level : %d\n",detailLevel);
1241 return TRUE;
1242 } else if ((key == 'L') || (key == 'l')) {
1243 detailLevel--;
1244 if (detailLevel < 0) detailLevel = 0;
1245 printf("level : %d\n",detailLevel);
1246 return TRUE;
1247 } else if ((key == 'b') || (key == 'B')) {
1248 drawType = 0;
1249 return TRUE;
1250 } else if ((key == 'd') || (key == 'D')) {
1251 drawType = 1;
1252 return TRUE;
1253 } else if ((key == 'p') || (key == 'P')) {
1254 drawType = 2;
1255 return TRUE;
1256 } else if ((key == 'q') || (key == 'Q')) {
1257 drawChannel--;
1258 if (drawChannel < 0) drawChannel = 0;
1259 printf("channel : %s\n",channels[drawChannel].name);
1260 return TRUE;
1261 } else if ((key == 'w') || (key == 'W')) {
1262 drawChannel++;
1263 if (drawChannel >= numChannels) drawChannel = numChannels-1;
1264 printf("channel : %s\n",channels[drawChannel].name);
1265 return TRUE;
1266 }
1267
1268 return FALSE;
1269 }
1270
1271
1272 ///////////////////////////////////////////////////////////////////////
1273 // Class : CBrickMap
1274 // Method : bound
1275 // Description : Bound the brickmap
1276 // Return Value : -
1277 // Comments :
1278 void CBrickMap::bound(float *bmin,float *bmax) {
1279 movvv(bmin,this->bmin);
1280 movvv(bmax,this->bmax);
1281 }
1282
1283
1284 ///////////////////////////////////////////////////////////////////////
1285 // Class : CBrickMap
1286 // Method : initBrickMap
1287 // Description : Initialize the brickmaps at the beginning of a frame
1288 // Return Value : -
1289 // Comments :
1290 void CBrickMap::initBrickMap(int m) {
1291 // This function is guaranteed to be called once and only once for each frame
1292 brickMaps = NULL;
1293 referenceNumber = 0;
1294 currentMemory = 0;
1295 maxMemory = m;
1296 }
1297
1298
1299
1300
1301 ///////////////////////////////////////////////////////////////////////
1302 // Class : CBrickMap
1303 // Method : flushBrickMap
1304 // Description : This function is called to re-claim some of the memory from the brickmap
1305 // Return Value : -
1306 // Comments :
1307 void CBrickMap::flushBrickMap(int allBricks) {
1308 int numNodes;
1309 CBrickNode **nodes;
1310 CBrickNode *cNode;
1311 CBrickMap *cMap;
1312 int i;
1313
1314 // Collect the loaded bricks into an array
1315 numNodes = 0;
1316 for (cMap=brickMaps;cMap!=NULL;cMap=cMap->nextMap) {
1317 for (i=0;i<BRICK_HASHSIZE;i++) {
1318 for (cNode=cMap->activeBricks[i];cNode!=NULL;cNode=cNode->next) {
1319 if (cNode->brick != NULL) numNodes++;
1320 }
1321 }
1322 }
1323
1324 nodes = new CBrickNode*[numNodes*2];
1325 numNodes = 0;
1326 for (cMap=brickMaps;cMap!=NULL;cMap=cMap->nextMap) {
1327 for (i=0;i<BRICK_HASHSIZE;i++) {
1328 for (cNode=cMap->activeBricks[i];cNode!=NULL;cNode=cNode->next) {
1329 if (cNode->brick != NULL) {
1330 nodes[numNodes*2] = cNode;
1331 nodes[numNodes*2+1] = (CBrickNode *) cMap;
1332 numNodes++;
1333 }
1334 }
1335 }
1336 }
1337
1338 // Sort the bricks wrt. to the last reference
1339 brickQuickSort(nodes,0,numNodes-1);
1340
1341 // Swap out the bricks
1342 if (allBricks == FALSE) {
1343 numNodes = numNodes >> 1;
1344 stats.numBrickmapCachePageouts += numNodes;
1345 }
1346
1347
1348 // Eliminate nodes
1349 for (i=0;i<numNodes;i++) {
1350 CBrickNode *cNode = nodes[i*2];
1351 CBrickMap *cMap = (CBrickMap *) nodes[i*2+1];
1352 CVoxel *cVoxel,*tVoxel;
1353 int j;
1354
1355 // Strategy - if we are modifying, save the contents to backing store
1356 // otherwise, just forget the voxel altogether
1357 if (cMap->modifying == TRUE) {
1358 // Write and free the brick
1359
1360 if (cNode->fileIndex == -1) {
1361 // If this is the first time we're writing, append it to the end
1362 fseek(cMap->file,0,SEEK_END);
1363 cNode->fileIndex = ftell(cMap->file);
1364 } else {
1365 // Go to the correct position
1366 fseek(cMap->file,cNode->fileIndex,SEEK_SET);
1367 }
1368
1369 uint32_t bs[BRICK_PRESENCE_LONGS];
1370
1371 // Write all voxels, do not spend time compressing
1372 for(j=0;j<BRICK_PRESENCE_LONGS;j++) bs[j] = 0xFFFFFFFFL;
1373
1374 fwrite(bs,sizeof(uint32_t)*BRICK_PRESENCE_LONGS,1,cMap->file);
1375
1376 for(j=BRICK_SIZE*BRICK_SIZE*BRICK_SIZE,cVoxel=cNode->brick->voxels;j>0;j--) {
1377 float *vdata = (float*) (cVoxel + 1);
1378
1379 fwrite(cVoxel,sizeof(CVoxel) + cMap->dataSize*sizeof(float),1,cMap->file);
1380
1381 while((tVoxel=cVoxel->next) != NULL) {
1382 cVoxel->next = tVoxel->next;
1383 fwrite(tVoxel,1,sizeof(CVoxel) + cMap->dataSize*sizeof(float),cMap->file);
1384 delete [] (char *) tVoxel;
1385
1386 currentMemory -= sizeof(CVoxel) + cMap->dataSize*sizeof(float);
1387 }
1388
1389 cVoxel = (CVoxel*) (vdata + cMap->dataSize);
1390 }
1391
1392 // Free the brick
1393 delete[] (char*) cNode->brick;
1394 cNode->brick = NULL;
1395
1396 // Update the used memory
1397 currentMemory -= sizeof(CBrick) + (sizeof(CVoxel) + cMap->dataSize*sizeof(float))*(BRICK_SIZE*BRICK_SIZE*BRICK_SIZE);
1398 } else {
1399 // Just free the brick
1400
1401 for(j=BRICK_SIZE*BRICK_SIZE*BRICK_SIZE,cVoxel=cNode->brick->voxels;j>0;j--) {
1402 float *vdata = (float*) (cVoxel + 1);
1403
1404 while((tVoxel=cVoxel->next) != NULL) {
1405 cVoxel->next = tVoxel->next;
1406 delete [] (char *) tVoxel;
1407
1408 currentMemory -= sizeof(CVoxel) + cMap->dataSize*sizeof(float);
1409 }
1410
1411 cVoxel = (CVoxel*) (vdata + cMap->dataSize);
1412 }
1413
1414 // Free the brick
1415 delete[] (char*) cNode->brick;
1416 cNode->brick = NULL;
1417
1418 // Update the used memory
1419 currentMemory -= sizeof(CBrick) + (sizeof(CVoxel) + cMap->dataSize*sizeof(float))*(BRICK_SIZE*BRICK_SIZE*BRICK_SIZE);
1420 }
1421 }
1422
1423 delete [] nodes;
1424 }
1425
1426
1427 ///////////////////////////////////////////////////////////////////////
1428 // Class : CBrickMap
1429 // Method : shutdownBrickMap
1430 // Description : Clear the memory used by the brickmaps
1431 // Return Value : -
1432 // Comments :
1433 void CBrickMap::shutdownBrickMap() {
1434 // This function is guaranteed to be called once and only once for each frame
1435 assert(currentMemory == 0);
1436 }
1437
1438
1439 ///////////////////////////////////////////////////////////////////////
1440 // Class : CBrickMap
1441 // Method : brickQuickSort
1442 // Description : Quick sort the bricks wrt. to the referenceNumbers
1443 // Return Value : -
1444 // Comments :
1445 void CBrickMap::brickQuickSort(CBrickNode **nodes,int start,int end) {
1446 int i,last;
1447 CBrickNode *cNode;
1448
1449 for (last=start,i=start+1;i<=end;i++) {
1450 if (nodes[i*2]->brick->referenceNumber < nodes[start*2]->brick->referenceNumber) {
1451 last++;
1452 cNode = nodes[last*2];
1453 nodes[last*2] = nodes[i*2];
1454 nodes[i*2] = cNode;
1455
1456 cNode = nodes[last*2+1];
1457 nodes[last*2+1] = nodes[i*2+1];
1458 nodes[i*2+1] = cNode;
1459 }
1460 }
1461
1462 cNode = nodes[last*2];
1463 nodes[last*2] = nodes[start*2];
1464 nodes[start*2] = cNode;
1465
1466 cNode = nodes[last*2+1];
1467 nodes[last*2+1] = nodes[start*2+1];
1468 nodes[start*2+1] = cNode;
1469
1470
1471
1472 // Speed is not an issue since this is not done very frequently, so recursion is OK
1473 if ((last-1) > start)
1474 brickQuickSort(nodes,start,last-1);
1475
1476 if (end > (last+1))
1477 brickQuickSort(nodes,last+1,end);
1478 }
1479
1480
1481
1482
1483
1484
1485
1486
1487
1488
1489
1490
1491
1492 ///////////////////////////////////////////////////////////////////////
1493 // Function : makeBrickMap
1494 // Description : This function creates the 3D baed texture from point cloud representation
1495 // Return Value : -
1496 // Comments :
1497 void makeBrickMap(int nb,const char **src,const char *dest,TSearchpath *searchPath,int n,const char **tokens,const void **params) {
1498 char tempName[OS_MAX_PATH_LENGTH];
1499 char fileName[OS_MAX_PATH_LENGTH];
1500 int i;
1501
1502 float maxVariation = 0.002f;
1503 float radiusScale = 1.0f;
1504 int maxDepth = 10;
1505 for(i =0;i<n;i++){
1506 if(!strcmp(tokens[i],"maxerror")){
1507 maxVariation = ((float*)params[i])[0];
1508 } else if (!strcmp(tokens[i],"radiusscale")){
1509 radiusScale = ((float*)params[i])[0];
1510 } else if (!strcmp(tokens[i],"maxdepth")){
1511 maxDepth = ((int*)params[i])[0];
1512 }
1513 }
1514
1515 // If not initialized already, init the brick memory manager
1516 // Use a large memory limit when creating brickmaps
1517 // Note: needed as RiMakeXYZ can only be called outside a frame, and then
1518 // the shading context is gone
1519 CBrickMap::initBrickMap(300000000);
1520
1521 // FIXME: deal with multiple brickmaps
1522 if (CRenderer::locateFile(fileName,src[0],searchPath)) {
1523 FILE *in;
1524 if ((in = ropen(fileName,"rb",filePointCloud,TRUE)) != NULL) {
1525
1526 // create backing store in a temp file
1527 // FIXME: make osTempname not always prefix dir
1528 sprintf(tempName,"%s.tmp",dest);
1529
1530 CPointCloud *cPtCloud = new CPointCloud(filePointCloud,identityMatrix,identityMatrix,in);
1531 CBrickMap *cBMap = new CBrickMap(tempName,cPtCloud->bmin,cPtCloud->bmax,identityMatrix,identityMatrix,cPtCloud->toNDC,cPtCloud->channels,cPtCloud->numChannels,maxDepth);
1532 float *data = cPtCloud->data.array;
1533 for (i=1;i<=cPtCloud->numItems;i++) {
1534 CPointCloudPoint *p = cPtCloud->items + i;
1535 float *C = data + p->entryNumber;
1536 const float R = p->dP * radiusScale;
1537 // guard against duff data making the map too deep
1538 if (R<C_EPSILON || R!=R) continue;
1539
1540 assert(inBox(cPtCloud->bmin,cPtCloud->bmax,p->P));
1541
1542 cBMap->store(C,p->P,p->N,R);
1543 }
1544
1545 cBMap->finalize();
1546
1547 // compact to the final file
1548 cBMap->compact(dest,maxVariation);
1549
1550 delete cBMap;
1551 delete cPtCloud;
1552 // clean up
1553 osDeleteFile(tempName);
1554 } else {
1555 error(CODE_BADTOKEN,"Point cloud file \"%s\" could not be opened\n");
1556 }
1557 } else {
1558 error(CODE_BADTOKEN,"Point cloud file \"%s\" not found\n");
1559 }
1560
1561 CBrickMap::shutdownBrickMap();
1562 }
1563
1564