~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

Linux Cross Reference
Tina5/tina-tools/tinatool/tlvision/tlvisTrn_procs.c

Version: ~
Architecture: ~ [ i386 ] ~ [ alpha ] ~ [ m68k ] ~ [ mips ] ~ [ ppc ] ~ [ sparc ] ~ [ sparc64 ] ~

  1 /**********
  2 *
  3 *
  4  * Copyright (c) 2003, Division of Imaging Science and Biomedical Engineering,
  5  * University of Manchester, UK.  All rights reserved.
  6  * 
  7  * Redistribution and use in source and binary forms, with or without modification, 
  8  * are permitted provided that the following conditions are met:
  9  * 
 10  *   . Redistributions of source code must retain the above copyright notice, 
 11  *     this list of conditions and the following disclaimer.
 12  *    
 13  *   . Redistributions in binary form must reproduce the above copyright notice,
 14  *     this list of conditions and the following disclaimer in the documentation 
 15  *     and/or other materials provided with the distribution.
 16  * 
 17  *   . Neither the name of the University of Manchester nor the names of its
 18  *     contributors may be used to endorse or promote products derived from this 
 19  *     software without specific prior written permission.
 20  * 
 21  * 
 22  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 
 23  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 
 24  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 
 25  * ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE 
 26  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 
 27  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 
 28  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
 29  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 
 30  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 
 31  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 
 32  * POSSIBILITY OF SUCH DAMAGE.
 33 **********
 34 *
 35 * Program :   TINA
 36 * File    :  $Source: /home/tina/cvs/tina-tools/tinatool/tlvision/tlvisTrn_procs.c,v $
 37 * Date    :  $Date: 2009/03/23 16:08:53 $
 38 * Version :  $Revision: 1.6 $
 39 * CVS Id  :  $Id: tlvisTrn_procs.c,v 1.6 2009/03/23 16:08:53 paul Exp $
 40 *
 41 * Notes   :
 42 *
 43 *
 44 *
 45 *********
 46 */
 47 
 48 
 49 #if HAVE_CONFIG_H
 50   #include <config.h>
 51 #endif
 52 
 53 #include <math.h>
 54 #include <stdlib.h>
 55 
 56 #include <tina/sys/sysDef.h>
 57 #include <tina/sys/sysPro.h>
 58 #include <tina/math/mathDef.h>
 59 #include <tina/math/mathPro.h>
 60 #include <tina/geometry/geomDef.h>
 61 #include <tina/geometry/geomPro.h>
 62 #include <tina/image/imgDef.h>
 63 #include <tina/image/imgPro.h>
 64 #include <tinatool/draw/drawPro.h>
 65 
 66 /* needed for triangle */
 67 #define REAL double
 68 
 69 #include <tinatool/tlvision/tlvisTrn_triangle.h>
 70 #include "tlvisTrn_procs.h"
 71 
 72 
 73 #define IMAGE_FILL 0
 74 #define SHADE_FILL 1
 75 #define MEAN_FILL  2
 76 
 77 
 78 static Terrain_data *terrain = NULL;
 79 static Imrect *im = NULL;
 80 static Imrect *shade_im = NULL;
 81 static Imrect *mean_im = NULL;
 82 
 83 static float scalez = 0.25;
 84 static int samplex = 16, sampley = 16;
 85 
 86 static Bool show_top = true;
 87 static Bool show_bot = false;
 88 static Bool show_lines = true;
 89 static Bool show_fill = false;
 90 static Bool show_all = false;
 91 
 92 static int fill_type = IMAGE_FILL;
 93 /* static int line_col; */
 94 
 95 void    terr_fill_set(int new_fill_type)
 96 {
 97     fill_type = new_fill_type;
 98 }
 99 
100 /*
101 void    terr_line_col_set(int col)
102 {
103     line_col = col;
104 }
105 */
106 
107 void    terr_show_lines_set(Bool show)
108 {
109     show_lines = show;
110 }
111 
112 void    terr_show_all_set(Bool show)
113 {
114     show_all = show;
115 }
116 
117 void    terr_show_fill_set(Bool show)
118 {
119     show_fill = show;
120 }
121 
122 void     terr_show(Bool top, Bool bot)
123 {
124     show_top = top;
125     show_bot = bot;
126 }
127 
128 void    terr_im_set(Imrect * new_im)
129 {
130     im = new_im;
131 }
132 
133 Imrect *terr_im_get(void)
134 {
135     return (im);
136 }
137 
138 void    terr_samplex_set(int sx)
139 {
140     samplex = sx;
141 }
142 
143 int     terr_samplex_get(void)
144 {
145     return (samplex);
146 }
147 
148 void    terr_sampley_set(int sy)
149 {
150     sampley = sy;
151 }
152 
153 int     terr_sampley_get(void)
154 {
155     return (sampley);
156 }
157 
158 void    terr_scalez_set(double sz)
159 {
160     scalez = sz;
161 }
162 
163 double  terr_scalez_get(void)
164 {
165     return (scalez);
166 }
167 
168 void    terrain_skeldraw(Tv * tv)
169 {
170     if (terrain == NULL)
171         return;
172 
173     tv_screen_double_buffer_start(tv->tv_screen);
174     tv_terrain_skel(tv, terrain->data, terrain->mask, terrain->m, terrain->n);
175     tv_screen_double_buffer_end(tv->tv_screen);
176 }
177 
178 static int patch_shade(int i, int j, Tv * tv)
179 {
180     Vec3    p00 = {Vec3_id};
181     Vec3    p01 = {Vec3_id};
182     Vec3    p10 = {Vec3_id};
183     Vec3    n = {Vec3_id};
184 
185     if (tv == NULL)
186         return (0);
187     p00 = terrain->data[i][j];
188     p01 = terrain->data[i][j + 1];
189     p10 = terrain->data[i + 1][j];
190     /* BUGFIX Vec3 p11; p11 = terrain->data[i + 1][j + 1]; */
191 
192     n = vec3_unit(vec3_cross(vec3_diff(p10, p00), vec3_diff(p01, p00)));
193     return ((int) (64.0 + 191.0 * fabs(vec3_dot(n, tv->ez3))));
194 }
195 
196 static int patch_mean(int i, int j)
197 {
198     double  sum = 0;
199 
200     sum += vec3_z(terrain->data[i][j]);
201     sum += vec3_z(terrain->data[i][j + 1]);
202     sum += vec3_z(terrain->data[i + 1][j]);
203     sum += vec3_z(terrain->data[i + 1][j + 1]);
204     sum /= 4.0 * scalez;
205     return ((sum > 255) ? 255 : sum);
206 }
207 
208 static Imrect *get_image_for_surf(int (*intensity_func) ( /* ??? */ ), void *data)
209 {
210     Imregion *region;
211     Imrect *im;
212     int     n, m;
213     int     i, j;
214     int     x, y;
215     int     lx, ux, ly, uy;
216     double  dx, dy;
217     int     g, *row;
218 
219     if (terrain == NULL)
220         return (NULL);
221 
222     m = terrain->m;
223     n = terrain->n;
224 
225     lx = vec3_x(terrain->data[0][0]);
226     ly = vec3_y(terrain->data[0][0]);
227     ux = vec3_x(terrain->data[m - 1][n - 1]) + 1;
228     uy = vec3_y(terrain->data[m - 1][n - 1]) + 1;
229     dx = (double) (ux - lx) / (n - 1);
230     dy = (double) (uy - ly) / (m - 1);
231 
232     region = roi_alloc(lx, ly, ux, uy);
233     im = im_alloc(uy, ux, region, uchar_v);
234     rfree((void *) region);
235 
236     row = ivector_alloc(lx, ux);
237     for (i = 0; i < m - 1; i++)
238     {
239         for (j = 0; j < n - 1; j++)
240         {
241             g = intensity_func(i, j, data);
242             for (x = j * dx; x < (j + 1) * dx && x + lx < ux; x++)
243                 row[lx + x] = g;
244         }
245         for (y = i * dy; y < (i + 1) * dy && y + ly < uy; y++)
246             im_put_row(row, im, y + ly, lx, ux);
247     }
248     ivector_free(row, lx);
249     return (im);
250 }
251 
252 void    terrain_fulldraw(Tv * tv)
253 {
254     int     n, m;
255 
256     if (terrain == NULL)
257         return;
258 
259     n = terrain->n;
260     m = terrain->m;
261 
262     /* tv_set_color(tv, line_col); */
263     (void) tv_set_op(tv, (show_lines == true) ? rop_copy : rop_null);
264     tv_terrain_hidden_init(tv, show_bot, show_top);
265     if (show_fill == true)
266     {
267         switch (fill_type)
268         {
269         case IMAGE_FILL:
270             tv_terrain_fill3(tv, terrain->data, m, n, im);
271             break;
272         case SHADE_FILL:
273             tv_terrain_fill3(tv, terrain->data, m, n, shade_im);
274             break;
275         case MEAN_FILL:
276             tv_terrain_fill3(tv, terrain->data, m, n, mean_im);
277             break;
278         }
279     } else if (show_all)
280         tv_terrain3(tv, terrain->data, m, n);
281     else
282         tv_terrain_hidden3(tv, terrain->data, m, n);
283     (void) tv_set_op(tv, rop_copy);
284 }
285 
286 void    terrain_backdraw(Tv * tv)
287 {
288 }
289 
290 void    terr_terrain_set(Terrain_data * new_terrain)
291 {
292     terrain_data_free(terrain);
293     terrain = new_terrain;
294     im_free(shade_im);
295     shade_im = get_image_for_surf(patch_shade, (void *) terrain_tv());
296     im_free(mean_im);
297     mean_im = get_image_for_surf(patch_mean, NULL);
298 }
299 
300 void    terrain_from_image(void)
301 {
302     if (im == NULL)
303         return;
304 
305     terrain_data_free(terrain);
306     terrain = im_surface(im, (Imrect *) NULL, samplex, sampley, scalez);
307     im_free(shade_im);
308     shade_im = get_image_for_surf(patch_shade, (void *) terrain_tv());
309     im_free(mean_im);
310     mean_im = get_image_for_surf(patch_mean, NULL);
311 }
312 
313 void    terrain_init(Tv * tv)
314 {
315     Vec3    centre = {Vec3_id};
316     Vec3    aim = {Vec3_id};
317     Vec3    down = {Vec3_id};
318     float   radius, pscale;
319     int     m, n;
320 
321     if (tv == NULL || terrain == NULL)
322         return;
323 
324     m = terrain->m;
325     n = terrain->n;
326     centre = terrain->data[m / 2][n / 2];
327     radius = 0.5 * vec3_mod(vec3_diff(terrain->data[0][0], terrain->data[m - 1][n - 1]));
328     pscale = 3;
329     /**
330     aim = vec3_unit(vec3(1.0, 1.0, -1.0));
331     down = vec3_minus(vec3_ey());
332     tv_set_axis(tv, vec3_minus(vec3_ez()));
333     **/
334     aim = vec3_unit(vec3(0.0, -0.5, 1.0));
335     down = vec3_ey();
336     tv_camera3(tv, centre, radius, pscale, aim, down);
337     tv_set_axis(tv, vec3_ez());
338 }
339 
340 Tv     *terrain_tv(void)
341 {
342     static Tv *tv = NULL;
343 
344     if (tv != NULL)
345         return (tv);
346 
347     tv = tv_create("Terrain");
348     (void) tv_set_backdraw(tv, terrain_backdraw);
349     (void) tv_set_fulldraw(tv, terrain_fulldraw);
350     (void) tv_set_skeldraw(tv, terrain_skeldraw);
351     tv_set_init(tv, terrain_init);
352     (void) tv_set_zoomlevel(tv, ZOOMGR);
353     return (tv);
354 }
355 
356 static void triangle_init(struct triangulateio *in, struct triangulateio *out)
357 {
358   /* initialise the input triangle structure */
359 
360   in->pointlist = NULL;
361   in->pointattributelist = NULL;
362   in->pointmarkerlist = NULL;
363   in->numberofpoints = 0;
364   in->numberofpointattributes = 0;
365 
366   in->trianglelist = NULL;
367   in->triangleattributelist = NULL;
368   in->trianglearealist = NULL;
369   in->neighborlist = NULL;
370   in->numberoftriangles = 0;
371   in->numberofcorners = 0;
372   in->numberoftriangleattributes = 0;
373 
374   in->segmentlist = NULL;
375   in->segmentmarkerlist = NULL;
376   in->numberofsegments = 0;
377 
378   in->holelist = NULL;
379   in->numberofholes = 0;
380 
381   in->regionlist = NULL;
382   in->numberofregions = 0;
383 
384   in->edgelist = NULL;
385   in->edgemarkerlist = NULL;
386   in->normlist = NULL;
387   in->numberofedges = 0;
388 
389   /* initialise the output triangle structure */
390 
391   out->pointlist = NULL;
392   out->pointattributelist = NULL;
393   out->pointmarkerlist = NULL;
394   out->numberofpoints = 0;
395   out->numberofpointattributes = 0;
396 
397   out->trianglelist = NULL;
398   out->triangleattributelist = NULL;
399   out->trianglearealist = NULL;
400   out->neighborlist = NULL;
401   out->numberoftriangles = 0;
402   out->numberofcorners = 0;
403   out->numberoftriangleattributes = 0;
404 
405   out->segmentlist = NULL;
406   out->segmentmarkerlist = NULL;
407   out->numberofsegments = 0;
408 
409   out->holelist = NULL;
410   out->numberofholes = 0;
411 
412   out->regionlist = NULL;
413   out->numberofregions = 0;
414 
415   out->edgelist = NULL;
416   out->edgemarkerlist = NULL;
417   out->normlist = NULL;
418   out->numberofedges = 0;
419 }
420 
421 /* same_side: returns true if pixel in question is on the
422    same side of the line (defined by p1, p2) as p3 */
423 
424 static Bool same_side(Vec2 p1, Vec2 p2, Vec2 p3, Vec2 pix)
425 {
426   float m, p3test, pixtest;
427   
428   /* vertical line case */
429   if (vec2_x(p1) == vec2_x(p2))
430     {
431       /* both p3 and pix to right of line */
432       if ((vec2_x(p3) > vec2_x(p1)) && (vec2_x(pix) >= vec2_x(p1)))
433         return(true);
434       /* both p3 and pix to left of line */
435       if ((vec2_x(p3) < vec2_x(p1)) && (vec2_x(pix) <= vec2_x(p1)))
436         return(true);
437     }
438 
439   /* horizontal line case */
440   else if (vec2_y(p1) == vec2_y(p2))
441     {
442       /* both p3 and pix above line */
443       if ((vec2_y(p3) > vec2_y(p1)) && (vec2_y(pix) >= vec2_y(p1)))
444         return(true);
445       /* both p3 and pix below line */
446       if ((vec2_y(p3) < vec2_y(p1)) && (vec2_y(pix) <= vec2_y(p1)))
447         return(true);
448     }
449 
450   /* sloping line case */
451   else
452     {
453       /* find the line gradient */
454       m = (vec2_y(p2) - vec2_y(p1)) / (vec2_x(p2) - vec2_x(p1));
455       
456       if (abs(m) <= 1)
457         {
458           /* use x knowns, look at y vals */
459           p3test = m * (vec2_x(p3) - vec2_x(p1)) + vec2_y(p1);
460           pixtest = m * (vec2_x(pix) - vec2_x(p1)) + vec2_y(p1);
461 
462           /* both points below line */
463           if ((p3test > vec2_y(p3)) && (pixtest >= vec2_y(pix)))
464             return(true);
465           /* both points above line */
466           if ((p3test < vec2_y(p3)) && (pixtest <= vec2_y(pix)))
467             return(true);
468         }
469       else
470         {
471           /* use y knowns, look at x vals */
472           p3test = vec2_x(p1) + (vec2_y(p3) - vec2_y(p1)) / m;
473           pixtest = vec2_x(p1) + (vec2_y(pix) - vec2_y(p1)) / m;
474 
475           /* both points to the left of line */
476           if ((p3test > vec2_x(p3)) && (pixtest >= vec2_x(pix)))
477             return(true);
478           /* both points to the right of line */
479           if ((p3test < vec2_x(p3)) && (pixtest <= vec2_x(pix)))
480             return(true);
481         }
482     }
483   
484   return(false);
485 }
486 
487 /* interp_val: given a triangle in 3d space (p1, p2, p3)
488    what is the z value of a point on that plane with x y
489    coordinates given by pix */
490 
491 static float interp_val(Vec3 p1, Vec3 p2, Vec3 p3, Vec2 pix)
492 {
493   Vec3 normal, interp;
494 
495   normal = vec3_cross(vec3_diff(p1, p3), vec3_diff(p2, p3));
496   interp = vec3_inter_line_plane(vec3(vec2_x(pix), vec2_y(pix), 0.0),
497                                  vec3(0.0, 0.0, 1.0),
498                                  p1, normal);
499   return(vec3_z(interp));
500 }
501 
502 /* triangulate_IM: uses delaunay triangulation routine to
503    divide up a set of points into triangles and then interpolates
504    those triangles to fill in the unknown points in the triangles */
505 /*
506 void triangulate_IM(Imrect *im)
507 {
508   struct triangulateio in, out;
509   Vec2 *lut = NULL, p1, p2, p3;
510   Vec3 t1, t2, t3;
511   int num, x, y, xmin, xmax, ymin, ymax;
512 
513   if (im == NULL || im->vtype!=float_v)
514     {
515       error("No IM or IM not of type float", non_fatal);
516       return;
517     }
518 
519   triangle_init(&in, &out);
520 
521   FOR_IM(im->region, y, x)
522     if (IM_FLOAT(im, y, x))
523       in.numberofpoints++;
524 
525   lut = (Vec2 *)calloc(in.numberofpoints, sizeof(Vec2));
526 
527 
528   in.pointlist = (REAL *)malloc(in.numberofpoints * 2 * sizeof(REAL));
529 
530   num = 0;
531 
532   FOR_IM(im->region, y, x)
533     {
534       if (IM_FLOAT(im, y, x) == 0.0)
535         continue;
536       
537       vec2_x(lut[num]) = x;
538       vec2_y(lut[num]) = y;
539       
540       in.pointlist[2*num] = x;
541       in.pointlist[2*num+1] = y;
542       num++;
543     }
544   
545   triangulate("zPNQ", &in, &out, NULL);
546 
547   for (num=0; num<out.numberoftriangles; num++)
548     {
549       p1 = lut[out.trianglelist[3*num]];
550       p2 = lut[out.trianglelist[3*num+1]];
551       p3 = lut[out.trianglelist[3*num+2]];
552 
553       xmin = tina_int(MIN(MIN(vec2_x(p1), vec2_x(p2)), vec2_x(p3)));
554       ymin = tina_int(MIN(MIN(vec2_y(p1), vec2_y(p2)), vec2_y(p3)));
555       xmax = tina_int(MAX(MAX(vec2_x(p1), vec2_x(p2)), vec2_x(p3)));
556       ymax = tina_int(MAX(MAX(vec2_y(p1), vec2_y(p2)), vec2_y(p3)));
557 
558       for (y=ymin; y<=ymax; y++)
559         for (x=xmin; x<=xmax; x++)
560           if ((IM_FLOAT(im, y, x) == 0.0) &&
561               same_side(p1, p2, p3, vec2(x, y)) &&
562               same_side(p1, p3, p2, vec2(x, y)) &&
563               same_side(p2, p3, p1, vec2(x, y)))
564             {
565               
566               t1 = vec3(vec2_x(p1), vec2_y(p1), 
567                         IM_FLOAT(im, tina_int(vec2_y(p1)), tina_int(vec2_x(p1))));
568               t2 = vec3(vec2_x(p2), vec2_y(p2), 
569                         IM_FLOAT(im, tina_int(vec2_y(p2)), tina_int(vec2_x(p2))));
570               t3 = vec3(vec2_x(p3), vec2_y(p3), 
571                         IM_FLOAT(im, tina_int(vec2_y(p3)), tina_int(vec2_x(p3))));
572               
573               IM_FLOAT(im, y, x) = interp_val(t1, t2, t3, vec2(x, y));
574             }
575     }
576   
577   cfree(lut);
578   free(in.pointlist);
579   free(out.trianglelist);
580 }
581 */
582 
583 /* triangulate_geom: uses delaunay triangulation routine to triangulate
584    a set of 3D points. the main triangulation routine is driven by
585    the x,y coordinates of the points. the results can then be displayed
586    in the threed tv for example */
587 
588 /* note: this proc adds the triangulation lines to the existing
589    geomlist, it does NOT free the given geomlist and return a
590    new one */
591 /*
592 List *triangulate_geom(List *geomlist)
593 {
594   struct triangulateio in, out;
595   List *ptr;
596   Point3 *point, **lut = NULL;
597   Vec3 p1, p2;
598   int num;
599 
600   if (geomlist == NULL)
601     return(geomlist);
602   
603   triangle_init(&in, &out);
604 
605   for (ptr=geomlist; ptr!=NULL; ptr=ptr->next)
606     if (ptr->type == POINT3)
607       in.numberofpoints++;
608 
609   lut = (Point3 **)calloc(in.numberofpoints, sizeof(Point3 *));
610 
611   in.pointlist = (REAL *)malloc(in.numberofpoints * 2 * sizeof(REAL));
612 
613   for (num=0, ptr=geomlist; ptr!=NULL; ptr=ptr->next)
614     {
615       if (ptr->type != POINT3)
616         continue;
617       
618       point = (Point3 *)ptr->to;
619 
620       lut[num] = point;
621       
622       in.pointlist[2*num] = vec3_x(point->p);
623       in.pointlist[2*num+1] = vec3_y(point->p);
624       num++;
625     }
626   
627   triangulate("zPNQ", &in, &out, NULL);
628 
629 
630   for (num=0; num<out.numberoftriangles; num++)
631     {
632       p1 = lut[out.trianglelist[3*num]]->p;
633       p2 = lut[out.trianglelist[3*num+1]]->p;
634       geomlist = list_addtostart(geomlist,
635                                  link_alloc((void *)line3_make(p1, p2, LINE3),
636                                             LINE3));
637 
638       p1 = lut[out.trianglelist[3*num+1]]->p;
639       p2 = lut[out.trianglelist[3*num+2]]->p;
640       geomlist = list_addtostart(geomlist,
641                                  link_alloc((void *)line3_make(p1, p2, LINE3),
642                                             LINE3));
643 
644       p1 = lut[out.trianglelist[3*num+2]]->p;
645       p2 = lut[out.trianglelist[3*num]]->p;
646       geomlist = list_addtostart(geomlist,
647                                  link_alloc((void *)line3_make(p1, p2, LINE3),
648                                             LINE3));
649     }
650   
651   cfree(lut);
652   free(in.pointlist);
653   free(out.trianglelist);
654 
655   return(geomlist);
656 }
657 */
658 

~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

This page was automatically generated by the LXR engine.
Visit the LXR main site for more information.