Initial population.

1 files changed, 675 insertions, 0 deletions

diff --git a/src/graph/view3D/contour3d.c.out b/src/graph/view3D/contour3d.c.out
new file mode 100755
index 00000000..14e791e2
--- /dev/null
+++ b/src/graph/view3D/contour3d.c.out
@@ -0,0 +1,675 @@
+  /***********************************************************************
+    contour.c
+    
+    begun 19 September 1992, Jim Wen
+   ***********************************************************************/
+
+#include "header.h"
+#include "draw.h"
+
+#define segmentDEBUG_X
+
+#define contourDEBUG_X
+#define use_min
+#define realloc_bug_fixed_NOT
+
+  /*=====================================================================*
+    Static variables
+   *=====================================================================*/
+int noo=0;
+poly   		       *contour_poly_list;
+poly     	       *active_list_first, *active_list_last, *active_list_current;
+segment_list_struct    *tmp_segment_list;
+
+  /*=====================================================================*
+    macro definitions
+   *=====================================================================*/
+#define foreach_poly(p)       for ((p)=contour_poly_list; p != NIL(poly); p=p->next)
+#define foreach_active(p,f,l) if (f != NIL(poly)) for (p=f; p!=l; p=(p)->next)
+#define in_range(x,a,b)	      ( ((x>=a) && (x<b)) || ((x<a) && (x>=b)) )
+#define poly_in_plane(p,z)    ( (p->contour_min==z) && (p->contour_max==z) )
+
+			       
+
+  /*=====================================================================*
+    local function declarations
+   *=====================================================================*/
+int 		 contour_compare();
+void 		 add_segment();
+segment_struct  *make_the_segment();
+viewTriple 	*mkpoint();
+CONTOUR_float    get_t_from_pts();
+void 		 make_active_list();
+int 		 maintain_active_list();
+void             contour_minMaxPolygons();
+
+  /*=====================================================================*
+    contour_compare()
+    
+    The compare function passed to msort. 
+   *=====================================================================*/
+int contour_compare(p1,p2)
+  poly *p1, *p2;
+{
+
+#ifdef use_min
+  if (p1->contour_min < p2->contour_min) return(-1);
+  else if (p1->contour_min == p2->contour_min) return(0);
+  else return(1);
+#else
+  if (p1->contour_max > p2->contour_max) return(-1);
+  else if (p1->contour_max == p2->contour_max) return(0);
+  else return(1);
+#endif /* use_min */
+}  /* contour_compare() */
+
+  /*=====================================================================*
+    do_contour_map()
+   *=====================================================================*/
+int do_contour_map()
+{
+
+  poly *pp;
+  poly *ap;
+  CONTOUR_float z_now;
+  int got_more;
+  int *anIndex;
+  viewTriple *aPt, *daPt, *one_point;
+  int jj, segment_index;
+  int got_one_intersection;
+  int done;
+  segment_list_struct *sl;
+#ifdef contour_object_rotate
+  float rotMat[4][4];
+  float transformed_zmax, transformed_zmin;
+#endif
+#ifdef conDEBUG
+  segment_struct *seg;
+#endif
+
+        /*---------------------------------------------------*
+	  the flag "did_contour" should be set to "no" whenever
+	  the user modifies one of the parameters
+         *---------------------------------------------------*/
+#ifdef contourDEBUG_x
+  fprintf(stderr,"Contour is %s\n",(did_contour)?"yes":"no");
+#endif
+  if (did_contour) return 1;
+  did_contour = yes;
+
+#ifdef contour_object_rotate
+        /*---------------------------------------------------*
+	  transform all the points for arbitrary plane
+	  slicing.  this includes all the viewTriples
+	  being referenced as well as the boundaries
+	  in viewData (to get contour_min, contour_max).
+         *---------------------------------------------------*/
+#ifdef oldie
+  rot_theta = 0;
+  rot_phi   = pi/2;
+#endif
+
+  sinTheta  = sin(-rot_theta);
+  cosTheta  = cos(-rot_theta);
+  sinPhi    = sin(rot_phi);
+  cosPhi    = cos(rot_phi);
+  ROTATE1(rotMat);
+
+        /*---------------------------------------------------*
+	  transform all the points.
+	  the zmin and zmax values need to be recalculated
+	  now that the object has been transformed.  note
+	  that transforming the z extreme points will not
+	  work correctly (i know, i've tried it).
+         *---------------------------------------------------*/
+    {
+      int i,j,k;
+      LLPoint *anLLPoint;
+      LPoint *anLPoint;
+      int *anIndex;
+      viewTriple *daPoint;
+      float v_in[4], v_out[4];
+      int first_time = yes;
+
+      anLLPoint = viewData.lllp.llp;
+      for (i=0; i<viewData.lllp.numOfComponents; i++,anLLPoint++) {
+	anLPoint = anLLPoint->lp;
+	for (j=0; j<anLLPoint->numOfLists; j++,anLPoint++) {
+	  anIndex = anLPoint->indices;
+	  for (k=0; k<anLPoint->numOfPoints; k++,anIndex++) {
+	    daPoint = refPt3D(viewData, *anIndex);
+        /*---------------------------------------------------*
+	  inefficient code so that the vector package could
+	  be used to see if things work;
+	  should change viewTriple's <x,y,z> to an array(?)
+         *---------------------------------------------------*/
+	    v_in[0] = daPoint->x;  v_in[1] = daPoint->y;  
+	    v_in[2] = daPoint->z;  v_in[3] = 1.0;
+	    vectorMatrix4(v_in, rotMat, v_out);
+	    daPoint->contour_x = v_out[0];  
+	    daPoint->contour_y = v_out[1];
+	    daPoint->contour_z = v_out[2];
+	    if (first_time) {
+	      first_time = no;
+	      transformed_zmin = transformed_zmax = v_out[2];
+	    } else {
+	      if (v_out[2] < transformed_zmin)
+	        transformed_zmin = v_out[2];
+	      else if (v_out[2] > transformed_zmax) 
+	        transformed_zmax = v_out[2];
+	    }
+	  } /* for points in LPoints (k) */
+	} /* for LPoints in LLPoints (j) */
+      } /* for LLPoints in LLLPoints (i) */
+
+    }
+
+#endif
+
+
+        /*---------------------------------------------------*
+	  set up the step size - it should be user adjustable
+	  and include all the slices possible.
+	  max_cuts is #slices seen
+	  z_step is #slices made
+	  cuts_used is #slices actually displayed
+         *---------------------------------------------------*/
+  z_step = (transformed_zmax - transformed_zmin)/max_cuts;
+#ifdef oldie
+  cuts_used = max_cuts;
+#endif
+
+        /*---------------------------------------------------*
+	  calculate the bounds of the polygons - the 
+	  contour_minMaxPolygons routine looks at the
+	  contour points rather than the object or 
+	  projected points
+         *---------------------------------------------------*/
+  contour_poly_list = copyPolygons(viewData.polygons);
+  contour_minMaxPolygons(contour_poly_list);
+
+        /*---------------------------------------------------*
+	  sort the polygons by the zmax value
+	  (or whatever transformed value, in general)
+         *---------------------------------------------------*/
+  contour_poly_list = msort(contour_poly_list, 0, 
+			    viewData.numPolygons, contour_compare);
+
+        /*---------------------------------------------------*
+	  having figured out how many cuts we need (should
+	  be a one time overhead so the following stuff
+	  should be in the initialization routine), we
+	  allocate an array of segment lists and initialize
+	  them (this part can stay here).
+         *---------------------------------------------------*/
+#ifdef oldie
+  segment_list = saymem("contour.c: segment_list", 
+			max_cuts, sizeof(segment_list_struct));
+#else
+        /*---------------------------------------------------*
+	  if the append flag is set, then we want to add 
+	  the new segment stuff onto the end of the old
+	  stuff (to build up a model of the surface).
+	  tmp_segment_list is use to keep track of the
+	  head of the new list while the routine goes
+	  through its paces.
+         *---------------------------------------------------*/
+  if (contour_append_lists_flag) {
+    noo = 0;
+#ifdef segmentDEBUG
+    fprintf(stderr,"======= series %d ========\n",noo);
+    fprintf(stderr,"     ---> before: sl->num=%d [%x]\n",
+	    segment_list->num, segment_list);
+#endif
+
+#ifdef realloc_bug_fixed
+    realloc(segment_list, (cuts_used + max_cuts) * sizeof(segment_list_struct));
+    tmp_segment_list = segment_list;
+    segment_list += cuts_used;  /* shift to end of list */
+    cuts_used += max_cuts;      /* size of new list */
+
+#ifdef segmentDEBUG
+    fprintf(stderr,"> %d cuts => seg at %x [old=%x]\n",
+  	  cuts_used, segment_list, tmp_segment_list);
+    fprintf(stderr,"        sl->num=%d, tsl->num=%d\n",
+	    segment_list->num, tmp_segment_list->num);
+#endif
+
+#else /* DONT_WORK */
+
+        /*---------------------------------------------------*
+	  Because realloc doesn't seem to work properly,
+	  we need to do this by hand
+         *---------------------------------------------------*/
+        /*---------------------------------------------------*
+	  allocate new space
+         *---------------------------------------------------*/
+  tmp_segment_list = saymem("contour.c: segment_list, append",
+			    cuts_used + max_cuts,
+			    sizeof(segment_list_struct));
+
+        /*---------------------------------------------------*
+	  copy over old data (1..cuts_used)
+         *---------------------------------------------------*/
+  {
+    segment_list_struct *tsl;
+
+    for (segment_index=0, sl=segment_list, tsl=tmp_segment_list; 
+	 segment_index<cuts_used; 
+	 segment_index++, sl++, tsl++) {
+      tsl->num = sl->num;
+      tsl->max_num = sl->max_num;
+      tsl->num_segs = sl->num_segs;
+      tsl->segments = sl->segments;
+    }
+  }
+
+        /*---------------------------------------------------*
+	  free the old stuff
+         *---------------------------------------------------*/
+/*  free(segment_list); */
+
+        /*---------------------------------------------------*
+	  now set segment_list to point to the point
+	  where tmp_segment_list stops - there ought
+	  to me max_cuts storage spaces left.
+         *---------------------------------------------------*/
+  segment_list = tmp_segment_list + cuts_used;
+
+        /*---------------------------------------------------*
+	  update cuts_used to have everything for a possible
+	  next iteration
+         *---------------------------------------------------*/
+  cuts_used += max_cuts;
+
+#endif /* realloc_BUG */
+
+
+  } else {
+    if (contour_allocated) {
+/*      free(segment_list); */
+    } else {
+      contour_allocated = yes;
+    }
+    noo = 0;
+    segment_list = saymem("contour.c: segment_list", 
+			  max_cuts, sizeof(segment_list_struct));
+    cuts_used = max_cuts;
+#ifdef segmentDEBUG
+    fprintf(stderr,"======= series %d ========\n",noo);
+    fprintf(stderr,"%d cuts => seg at %x\n",cuts_used, segment_list);
+#endif
+  }
+#endif
+
+  for (segment_index=0, sl=segment_list; 
+       segment_index<max_cuts; 
+       segment_index++, sl++) {
+    sl->num = ++noo;
+    sl->max_num = max_cuts;
+#ifdef segmentDEBUG
+    fprintf(stderr,"Made segment list %d  [%x]\n",noo,sl);
+/*    fprintf(stderr,"       ...(tsl->num=%d [%x]\n",
+	    tmp_segment_list->num, tmp_segment_list); */
+#endif
+    sl->num_segs = 0;
+    sl->segments = NIL(segment_struct);
+  }
+        /*---------------------------------------------------*
+	  create an "active_list" of polygons such that 
+	  this slice step intersects all and only those
+	  polygons in the active list.
+         *---------------------------------------------------*/
+  make_active_list(contour_poly_list,
+#ifdef use_min
+		   transformed_zmin, transformed_zmin+z_step,
+#else
+		   transformed_zmax, transformed_zmax-z_step,
+#endif
+		   &active_list_first, &active_list_last, &active_list_current);
+
+        /*---------------------------------------------------*
+	  iterate from zmax down to zmin with z_step increments
+         *---------------------------------------------------*/
+  segment_index = 0;
+#ifdef use_min
+  for (z_now=transformed_zmin; z_now<transformed_zmax; /* see below for incr*/) {
+#else
+  for (z_now=transformed_zmax; z_now>transformed_zmin; /* see below for incr*/) {
+#endif
+        /*---------------------------------------------------*
+	  for each of the polygons on the active list, 
+	  intersect each of line equation for the sides	
+	  with the plane equation for the plane at z_now.
+	  one of the following may occur:
+	     no intersections: haha - can't happen coz active
+	    one intersection : at the point, create point line
+	    two intersections: create line connecting points
+	    lies in the plane: create three segments
+	  note that this is a fairly inefficient approach but
+	  i'm just throwing this stuff together in an afternoon
+	  to see how it looks
+         *---------------------------------------------------*/
+    foreach_active(ap, active_list_first, active_list_last) {
+
+        /*---------------------------------------------------*
+	  do line-plane equation for each side of the 
+	  polygon
+	  for now - just 3+ sided polygons (no degenerates)
+         *---------------------------------------------------*/
+
+      if (ap->numpts >= 3) {
+	if (poly_in_plane(ap, z_now)) {
+        /*---------------------------------------------------*
+	  re-create all the segments of the polygon
+         *---------------------------------------------------*/
+	  daPt = refPt3D(viewData, *(ap->indexPtr + (ap->numpts - 1)));
+	  for (jj=0, anIndex=ap->indexPtr; jj<ap->numpts; jj++, anIndex++) {
+	    aPt = refPt3D(viewData, *anIndex);
+	    add_segment(segment_list, segment_index, 
+			make_the_segment(aPt, daPt));
+	    daPt = aPt;
+	  }
+	} else {
+        /*---------------------------------------------------*
+	  find the line that defines the intersection of 
+	  the polygon with the z-plane
+         *---------------------------------------------------*/
+	  got_one_intersection = no;
+	  done = no;
+	  daPt = refPt3D(viewData, *(ap->indexPtr + (ap->numpts - 1)));
+	  for (jj=0, anIndex=ap->indexPtr; 
+	       !done && jj<ap->numpts; 
+	       jj++, anIndex++) {
+	    aPt = refPt3D(viewData, *anIndex);
+	    if (in_range(z_now, aPt->contour_z, daPt->contour_z)) {
+	      if (got_one_intersection) {
+		add_segment(segment_list, segment_index,
+			    make_the_segment(one_point, 
+					     mkpoint(aPt, daPt, z_now)));
+		done = yes;
+	      } else {
+		one_point = mkpoint(aPt, daPt, z_now);
+		got_one_intersection = yes;
+	      }
+	    } 
+	    daPt = aPt;
+	  }  /* for */
+	}  /* else not lie in plane */
+
+      }
+    }  /* foreach_active(ap) */
+
+        /*---------------------------------------------------*
+	  maintain/update the active list, pruning off things
+	  the fall off the top (beyond z_now - z_step) and
+===>	  adding on things to the bottom that now fall into   
+	  the range [z_now ---> z_now-z_step].
+         *---------------------------------------------------*/
+    segment_index++;
+#ifdef use_min
+    z_now += z_step;
+    got_more = maintain_active_list(&active_list_first, 
+				    &active_list_last, 
+				    &active_list_current,
+				    z_now + z_step, 
+				    z_now);
+#else
+    z_now -= z_step;
+    got_more = maintain_active_list(&active_list_first, 
+				    &active_list_last, 
+				    &active_list_current,
+				    z_now, 
+				    z_now - z_step);
+#endif
+  }  /* for z_now from zmax to zmin */
+
+        /*---------------------------------------------------*
+	  if the segment lists have been appended, reset
+	  the global segment lists pointer to the top of
+	  the lists
+         *---------------------------------------------------*/
+  if (contour_append_lists_flag) {
+    segment_list = tmp_segment_list;
+#ifdef segmentDEBUG
+#ifdef oldie
+    fprintf(stderr,"    setting seg to old=%x\n", segment_list);
+    fprintf(stderr,"        num is %d\n",segment_list->num);
+#endif
+  {
+
+    for (segment_index=0, sl=segment_list; 
+	 segment_index<cuts_used; 
+	 segment_index++, sl++) {
+      fprintf(stderr,"   sl->num = %d\n",sl->num);
+    }
+  }
+
+#endif
+  }
+
+}  /* do_contour_map() */
+
+
+  /*=====================================================================*
+    make_active_list(da_list, z_min, z_max, af, al, ac)
+   *=====================================================================*/
+void make_active_list(da_list, z_min, z_max, af, al, ac)
+  poly *da_list;
+  CONTOUR_float z_min, z_max;
+  poly **af, **al, **ac;
+{
+
+  poly *tmp_p;
+
+        /*---------------------------------------------------*
+	  the first active polygon is the first one in the
+	  given, sorted list.  note that if it doesn't fall
+===>	  inside the z_max --> z_min range, af is set to NIL.
+         *---------------------------------------------------*/
+#ifdef use_min
+  if (da_list->contour_min > z_max) {
+#else
+  if (da_list->contour_max < z_min) {
+#endif
+    *af = NIL(poly);
+    return;
+  } else {
+    *af = da_list;
+  }
+
+        /*---------------------------------------------------*
+	  the current active polygon is set to "af" at this
+	  point but it could be the case that it is set to	
+	  "al" if, for example, the current span of z-values
+	  has no polygons - af is set to NIL, and the next
+	  time we look we start at ac=af.
+         *---------------------------------------------------*/
+  *ac = *af;
+
+        /*---------------------------------------------------*
+	  the last active polygon is the polygon right before
+===>	  the first one whose zmax is too small to make the
+	  list 
+         *---------------------------------------------------*/
+  *al = da_list;	
+  tmp_p = da_list->next;
+  for (; tmp_p != NIL(poly); tmp_p = tmp_p->next) {
+#ifdef use_min
+    if (tmp_p->contour_min > z_max) return;
+#else
+    if (tmp_p->contour_max < z_min) return;
+#endif
+    *al = tmp_p;
+  }
+  
+}  /* make_active_list() */
+
+
+  /*=====================================================================*
+    maintain_active_list(af, al, ac, z_max, z_min)
+   *=====================================================================*/
+int maintain_active_list(af, al, ac, z_max, z_min)
+  poly **al, **af, **ac;
+  CONTOUR_float z_max, z_min;
+{
+
+  poly *tmp_p;
+
+        /*---------------------------------------------------*
+	  first, get the lower boundary to be within range,
+	  pruning elements from the head of the list
+         *---------------------------------------------------*/
+  *af = *ac;
+#ifdef use_min
+  while ((*af) && (*af)->contour_max < z_min) {
+#else
+  while ((*af) && (*af)->contour_min > z_max) {
+#endif
+    *af = (*af)->next;
+  }
+
+        /*---------------------------------------------------*
+	  check to see if the upper boundary is still in
+	  range
+         *---------------------------------------------------*/
+#ifdef use_min
+  if ((*af) == NIL(poly) || ((*af)->contour_min > z_max)) {
+#else
+  if ((*af) == NIL(poly) || ((*af)->contour_max < z_min)) {
+#endif
+    /*   --- nope, it wasn't ---   */
+    *ac = *af;
+    *af = NIL(poly);
+    return(0);
+  }
+
+        /*---------------------------------------------------*
+	  upper boundary is okay...see if we need to add to
+	  the list on the lower bound side
+         *---------------------------------------------------*/
+  tmp_p = (*al)->next;
+  for (; tmp_p != NIL(poly); tmp_p = tmp_p->next) {
+#ifdef use_min
+    if (tmp_p->contour_min > z_max) return;
+#else
+    if (tmp_p->contour_max < z_min) return;
+#endif
+    *al = tmp_p;
+  }
+
+  return(1);
+
+}  /* maintain_active_list() */
+
+
+  /*=====================================================================*
+    add_segment(seg_list, index, seg)
+   *=====================================================================*/
+void add_segment(seg_list, index, seg)
+  segment_list_struct *seg_list;
+  int	index;
+  segment_struct *seg;
+{
+
+  segment_list_struct *sl;
+
+  sl = seg_list + index;
+
+  seg->next = sl->segments;
+  sl->segments = seg;
+  sl->num_segs++;
+
+}  /* add_segment() */
+
+
+  /*=====================================================================*
+    make_the_segment(pt1, pt2)
+   *=====================================================================*/
+segment_struct *make_the_segment(pt1, pt2)
+  viewTriple *pt1, *pt2;
+{
+  segment_struct *seg;
+  
+  seg = (segment_struct *)saymem("contour.c: segment",1,sizeof(segment_struct));
+  seg->point1 = pt1;
+  seg->point2 = pt2;
+
+  return(seg);
+
+}  /* make_the_segment() */
+
+
+  /*=====================================================================*
+    viewTriple *mkpoint(vt1, vt2, z_val)
+   *=====================================================================*/
+viewTriple *mkpoint(vt1, vt2, z_val)
+  viewTriple *vt1, *vt2;
+  CONTOUR_float z_val;
+{
+
+  viewTriple *vt;
+  CONTOUR_float t;
+
+  vt = (viewTriple *)saymem("contour.c: viewTriple",1,sizeof(viewTriple));
+  
+  t = get_t_from_pts(vt1->contour_z, vt2->contour_z, z_val);
+
+#ifdef waitaminute
+  vt->x = vt1->contour_x + (vt2->contour_x - vt1->contour_x) * t;
+  vt->y = vt1->contour_y + (vt2->contour_y - vt1->contour_y) * t;
+  vt->z = z_val;
+#else
+  vt->x = vt1->x + (vt2->x - vt1->x) * t;
+  vt->y = vt1->y + (vt2->y - vt1->y) * t;
+  vt->z = vt1->z + (vt2->z - vt1->z) * t;
+#endif
+
+  vt->contour_x = vt1->contour_x + (vt2->contour_x - vt1->contour_x) * t;
+  vt->contour_y = vt1->contour_y + (vt2->contour_y - vt1->contour_y) * t;
+  vt->contour_z = z_val;
+
+  return(vt);
+
+}  /* mkpoint() */
+
+
+
+  /*=====================================================================*
+    get_t_from_pts(z_min, z_max, z_val)
+   *=====================================================================*/
+CONTOUR_float get_t_from_pts(z_min, z_max, z_val)
+  CONTOUR_float z_min, z_max, z_val;
+{
+  CONTOUR_float t;
+
+  if (z_min == z_max) return 0;
+  
+  t = (z_val - z_min)/(z_max - z_min);
+
+  return(t);
+}  /* get_t_from_pts() */
+
+
+void contour_minMaxPolygons(aPoly)
+  poly *aPoly;
+{
+
+  int *anIndex;
+  int i;
+
+
+  for (; aPoly != NIL(poly); aPoly = aPoly->next) {
+    anIndex = aPoly->indexPtr;
+    aPoly->contour_min = aPoly->contour_max = 
+      refPt3D(viewData,*anIndex)->contour_z;
+    for (i=1,anIndex++; i<aPoly->numpts; i++,anIndex++) {
+      if (refPt3D(viewData,*anIndex)->contour_z < aPoly->contour_min)
+        aPoly->contour_min = refPt3D(viewData,*anIndex)->contour_z;
+      else if (refPt3D(viewData,*anIndex)->contour_z > aPoly->contour_max)
+        aPoly->contour_max = refPt3D(viewData,*anIndex)->contour_z;
+    }
+  }
+
+
+}  /* contour_minMaxPolygons */