ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/warp3d.c
(Generate patch)

Comparing ray/src/px/warp3d.c (file contents):
Revision 3.1 by greg, Tue Feb 4 16:03:48 1997 UTC vs.
Revision 3.5 by greg, Sat Feb 22 02:07:28 2003 UTC

# Line 1 | Line 1
1 /* Copyright (c) 1997 Regents of the University of California */
2
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ LBL";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6
4   /*
5   * 3D warping routines.
6   */
7  
8   #include <stdio.h>
9 <
9 > #include <stdlib.h>
10   #include <math.h>
14
11   #include "fvect.h"
16
12   #include "warp3d.h"
13  
14   #define MIND    1e-5            /* minimum distance between input points */
15  
16   typedef struct {
17          GNDX    n;              /* index must be first */
18 <        W3VEC   p;
18 >        W3VEC   v;
19   } KEYDP;                /* key/data allocation pair */
20  
21 < #define fgetvec(f,v)    (fgetval(f,'f',v) > 0 && fgetval(f,'f',v+1) > 0 \
22 <                                && fgetval(f,'f',v+2) > 0)
21 > #define fgetvec(p,v)    (fgetval(p,'f',v) > 0 && fgetval(p,'f',v+1) > 0 \
22 >                                && fgetval(p,'f',v+2) > 0)
23  
24 < #define AHUNK   32              /* number of points to allocate at a time */
24 > #define AHUNK   24              /* number of points to allocate at a time */
25  
31 #ifndef malloc
32 extern char     *malloc(), *realloc();
33 #endif
34 extern void     free();
26  
36
27   double
28   wpdist2(p1, p2)                 /* compute square of distance between points */
29   register W3VEC  p1, p2;
# Line 54 | Line 44 | register WARP3D        *wp;
44   {
45          int     rval = W3OK;
46          GNDX    gi;
47 <        W3VEC   gd;
47 >        W3VEC   gd, ov;
48  
49 <        if (wp->grid.gn[0] == 0 && (rval = new3dwgrid(wp)) != W3OK)
49 >        if (wp->grid.gn[0] == 0 && (rval = new3dgrid(wp)) != W3OK)
50                  return(rval);
51          rval |= gridpoint(gi, gd, pi, &wp->grid);
52 <        if (wp->grid.flags & W3EXACT) {
53 <                rval |= warp3dex(po, pi, wp);
54 <                return(rval);
55 <        }
56 <        if (wp->grid.flags & W3FAST) {
57 <                rval |= get3dgpt(po, gi, wp);
58 <                return(rval);
59 <        }
60 <        rval |= get3dgin(po, gi, gd, wp);
52 >        if (wp->grid.flags & W3EXACT)
53 >                rval |= warp3dex(ov, pi, wp);
54 >        else if (wp->grid.flags & W3FAST)
55 >                rval |= get3dgpt(ov, gi, wp);
56 >        else
57 >                rval |= get3dgin(ov, gi, gd, wp);
58 >        po[0] = pi[0] + ov[0];
59 >        po[1] = pi[1] + ov[1];
60 >        po[2] = pi[2] + ov[2];
61          return(rval);
62   }
63  
# Line 98 | Line 88 | register struct grid3d *gp;
88  
89  
90   int
91 < get3dgpt(po, ndx, wp)           /* get value for voxel */
92 < W3VEC   po;
91 > get3dgpt(ov, ndx, wp)           /* get value for voxel */
92 > W3VEC   ov;
93   GNDX    ndx;
94   register WARP3D *wp;
95   {
# Line 109 | Line 99 | register WARP3D        *wp;
99          int     rval = W3OK;
100          register int    i;
101  
102 <        le = lu_find(&wp->grid.gtab, (char *)ndx);
102 >        le = lu_find(&wp->grid.gtab, ndx);
103          if (le == NULL)
104                  return(W3ERROR);
105          if (le->data == NULL) {
# Line 123 | Line 113 | register WARP3D        *wp;
113                          if (wp->grid.flags & W3FAST)    /* on centers */
114                                  gpt[i] += .5*wp->grid.gstep[i];
115                  }
116 <                rval |= warp3dex(kd->p, gpt, wp);
116 >                rval = warp3dex(kd->v, gpt, wp);
117                  le->key = (char *)kd->n;
118 <                le->data = (char *)kd->p;
118 >                le->data = (char *)kd->v;
119          }
120 <        W3VCPY(po, (float *)le->data);
120 >        W3VCPY(ov, (float *)le->data);
121          return(rval);
122   }
123  
124  
125   int
126 < get3dgin(po, ndx, rem, wp)      /* interpolate from warp grid */
127 < W3VEC   po, rem;
126 > get3dgin(ov, ndx, rem, wp)      /* interpolate from warp grid */
127 > W3VEC   ov, rem;
128   GNDX    ndx;
129   WARP3D  *wp;
130   {
131          W3VEC   cv[8];
142        int     rval;
132          GNDX    gi;
133 +        int     rval = W3OK;
134          register int    i;
135                                          /* get corner values */
136          for (i = 0; i < 8; i++) {
# Line 151 | Line 141 | WARP3D *wp;
141          }
142          if (rval & W3ERROR)
143                  return(rval);
144 <        l3interp(po, cv, rem, 3);       /* perform trilinear interpolation */
144 >        l3interp(ov, cv, rem, 3);       /* perform trilinear interpolation */
145          return(rval);
146   }
147  
# Line 176 | Line 166 | int    n;
166  
167  
168   int
169 < warp3dex(po, pi, wp)            /* compute warp using 1/r^2 weighted avg. */
170 < W3VEC   po, pi;
169 > warp3dex(ov, pi, wp)            /* compute warp using 1/r^2 weighted avg. */
170 > W3VEC   ov, pi;
171   register WARP3D *wp;
172   {
173          double  d2, w, wsum;
174 <        W3VEC   pt;
174 >        W3VEC   vt;
175          register int    i;
176  
177 <        pt[0] = pt[1] = pt[2] = 0.;
177 >        vt[0] = vt[1] = vt[2] = 0.;
178          wsum = 0.;
179          for (i = wp->npts; i--; ) {
180                  d2 = wpdist2(pi, wp->ip[i]);
181                  if (d2 <= MIND*MIND) w = 1./(MIND*MIND);
182                  else w = 1./d2;
183 <                pt[0] += w*wp->op[i][0];
184 <                pt[1] += w*wp->op[i][1];
185 <                pt[2] += w*wp->op[i][2];
183 >                vt[0] += w*wp->ov[i][0];
184 >                vt[1] += w*wp->ov[i][1];
185 >                vt[2] += w*wp->ov[i][2];
186                  wsum += w;
187          }
188          if (wsum > 0.) {
189 <                po[0] = pt[0]/wsum;
190 <                po[1] = pt[1]/wsum;
191 <                po[2] = pt[2]/wsum;
189 >                ov[0] = vt[0]/wsum;
190 >                ov[1] = vt[1]/wsum;
191 >                ov[2] = vt[2]/wsum;
192          }
193          return(W3OK);
194   }
# Line 218 | Line 208 | int    flgs;
208                  eputs("new3dw: no memory\n");
209                  return(NULL);
210          }
211 <        wp->ip = wp->op = NULL;
211 >        wp->ip = wp->ov = NULL;
212          wp->npts = 0;
213          wp->grid.flags = flgs;
214          wp->grid.gn[0] = wp->grid.gn[1] = wp->grid.gn[2] = 0;
# Line 239 | Line 229 | WARP3D *wp;
229                  eputs(": cannot open\n");
230                  return(NULL);
231          }
232 <        if (wp == NULL && (wp = new3dw()) == NULL)
232 >        if (wp == NULL && (wp = new3dw(0)) == NULL)
233                  goto memerr;
234          while (fgetvec(fp, inp) && fgetvec(fp, outp))
235                  if (!add3dpt(wp, inp, outp))
# Line 258 | Line 248 | cleanup:
248  
249  
250   int
251 + set3dwfl(wp, flgs)              /* reset warp flags */
252 + register WARP3D *wp;
253 + int     flgs;
254 + {
255 +        if (flgs == wp->grid.flags)
256 +                return(0);
257 +        if ((flgs & (W3EXACT|W3FAST)) == (W3EXACT|W3FAST)) {
258 +                eputs("set3dwfl: only one of W3EXACT or W3FAST\n");
259 +                return(-1);
260 +        }
261 +        wp->grid.flags = flgs;
262 +        done3dgrid(&wp->grid);          /* old grid is invalid */
263 +        return(0);
264 + }
265 +
266 +
267 + int
268   add3dpt(wp, pti, pto)           /* add 3D point pair to warp structure */
269   register WARP3D *wp;
270   W3VEC   pti, pto;
# Line 269 | Line 276 | W3VEC  pti, pto;
276          if (wp->npts == 0) {                    /* initialize */
277                  wp->ip = (W3VEC *)malloc(AHUNK*sizeof(W3VEC));
278                  if (wp->ip == NULL) return(0);
279 <                wp->op = (W3VEC *)malloc(AHUNK*sizeof(W3VEC));
280 <                if (wp->op == NULL) return(0);
279 >                wp->ov = (W3VEC *)malloc(AHUNK*sizeof(W3VEC));
280 >                if (wp->ov == NULL) return(0);
281                  wp->d2min = 1e10;
282                  wp->d2max = 0.;
283                  W3VCPY(wp->llim, pti);
# Line 281 | Line 288 | W3VEC  pti, pto;
288                                          (wp->npts+AHUNK)*sizeof(W3VEC));
289                          if (na == NULL) return(0);
290                          wp->ip = na;
291 <                        na = (W3VEC *)realloc((char *)wp->op,
291 >                        na = (W3VEC *)realloc((char *)wp->ov,
292                                          (wp->npts+AHUNK)*sizeof(W3VEC));
293                          if (na == NULL) return(0);
294 <                        wp->op = na;
294 >                        wp->ov = na;
295                  }
296                  for (i = 0; i < 3; i++)         /* check boundaries */
297                          if (pti[i] < wp->llim[i])
# Line 302 | Line 309 | W3VEC  pti, pto;
309                  }
310          }
311          W3VCPY(wp->ip[wp->npts], pti);  /* add new point */
312 <        W3VCPY(wp->op[wp->npts], pto);
312 >        wp->ov[wp->npts][0] = pto[0] - pti[0];
313 >        wp->ov[wp->npts][1] = pto[1] - pti[1];
314 >        wp->ov[wp->npts][2] = pto[2] - pti[2];
315          done3dgrid(&wp->grid);          /* old grid is invalid */
316          return(++wp->npts);
317   }
# Line 312 | Line 321 | free3dw(wp)                    /* free WARP3D data */
321   register WARP3D *wp;
322   {
323          done3dgrid(&wp->grid);
324 <        free((char *)wp->ip);
325 <        free((char *)wp->op);
326 <        free((char *)wp);
324 >        free((void *)wp->ip);
325 >        free((void *)wp->ov);
326 >        free((void *)wp);
327   }
328  
329  
330 < long
330 > unsigned long
331   gridhash(gp)                    /* hash a grid point index */
332   GNDX    gp;
333   {
# Line 327 | Line 336 | GNDX   gp;
336  
337  
338   int
339 < new3dwgrid(wp)                  /* initialize interpolating grid for warp */
339 > new3dgrid(wp)                   /* initialize interpolating grid for warp */
340   register WARP3D *wp;
341   {
342 +        W3VEC   gmax;
343          double  gridstep, d;
344 +        int     n;
345          register int    i;
346                                  /* free old grid (if any) */
347          done3dgrid(&wp->grid);
# Line 341 | Line 352 | register WARP3D        *wp;
352                  return(W3BADMAP);       /* not enough disparity */
353                                  /* compute gamut */
354          W3VCPY(wp->grid.gmin, wp->llim);
355 <        W3VCPY(wp->grid.gmax, wp->ulim);
355 >        W3VCPY(gmax, wp->ulim);
356          gridstep = sqrt(wp->d2min);
357          for (i = 0; i < 3; i++) {
358                  wp->grid.gmin[i] -= .501*gridstep;
359 <                wp->grid.gmax[i] += .501*gridstep;
359 >                gmax[i] += .501*gridstep;
360          }
361          if (wp->grid.flags & W3EXACT) {
362                  wp->grid.gn[0] = wp->grid.gn[1] = wp->grid.gn[2] = 1;
363 +                wp->grid.gstep[0] = gmax[0] - wp->grid.gmin[0];
364 +                wp->grid.gstep[1] = gmax[1] - wp->grid.gmin[1];
365 +                wp->grid.gstep[2] = gmax[2] - wp->grid.gmin[2];
366                  return(W3OK);           /* no interpolation, so no grid */
367          }
368                                  /* create grid */
369          for (i = 0; i < 3; i++) {
370 <                d = wp->grid.gmax[i] - wp->grid.gmin[i];
371 <                wp->grid.gn[i] = d/gridstep + .5;
372 <                if (wp->grid.gn[i] >= MAXGN)
373 <                        wp->grid.gn[i] = MAXGN-1;
374 <                wp->grid.gstep[i] = d / wp->grid.gn[i];
370 >                d = gmax[i] - wp->grid.gmin[i];
371 >                n = d/gridstep + .5;
372 >                if (n >= MAXGN-1)
373 >                        n = MAXGN-2;
374 >                wp->grid.gstep[i] = d / n;
375 >                wp->grid.gn[i] = n;
376          }
377                                  /* initialize lookup table */
378          wp->grid.gtab.hashf = gridhash;

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines