#define PROG_NAME "SelectSubdivision"
#define PROG_DESC "???"
#define PROG_VERS "1.0"

/* Last edited on 2007-02-03 23:24:03 by stolfi */

#define SelectSubdivision_C_COPYRIGHT \
  ""
  
#define PROG_INFO \
  "" \
  " "


/*
  This module implements the "selective subdivision" of the k-elements
  ( 1 <= k <=3)  of a given 3D cellular map that are involved in geometric
  degeneracies.

  The option "element" allows to choose the selective subdivision of the
  k-element. The selective subdivision of a k-element also propagates the
  subdivision of their star (neighbors elements). */ 

#define SelectSubdivision_C_author \
  "Implemented by L. P. Lozada.\n" \
  "Revisions:\n" \
  "  08-06-2000:  Fixed and removed a bug that increases the number of\n" \
  "               original elements.\n" \
  "  11-06-2000:  Added the heredity of the root attributes for edges\n" \
  "               and walls, and we optimize the code. Hide the new elements\n" \
  "               insert in the refiment process:\n" \
  "               walls, edges and nodes are set with the attibute \n" \
  "               {exists=FALSE}."

#include <Mis.h>
#include <string.h>
#include <r3.h>
#include <Triangulation.h>
#define _GNU_SOURCE
#include <stdio.h>
#include <argparser.h>
#include <Octf.h>

uint NVE = 0;
  NVF: uint = 0;
  NVP: uint = 0;
                                 
TYPE
double Element = { @{Edge->?}, Wall, Tetrahedron };

typedef struct Options_t {
    char *inFile;            /* Initial guess file name (minus ".tp") */
    char *outFile;            /* Output file name prefix */
    Element element;
    char *elename;
  }

Options_t *GetOptions(int argc, char **argv);

int main(int argc, char **argv)
  
  { Options_t *o = GetOptions(argc, argv);
    char *topo_cmt = jsprintf("Created by %s on %s", PROG_NAME, Today());
    /* Random_t coins = MakeRandomSource(4615); */
    ??? tc = Triangulation.ReadToMa(o->inFile);
    ??? top = tc.top;
    {
      if (o->element == Element.edge)
      {
        /* Rescue the overall attributes of the triangulation */
        for (i = 0; i < top->wall.nelE; i++)
        {
          ??? p = top->wedge[i];
          Edge_t e = PEdge(p);
          with (en == e->num,
          double f = PWall(p), fn == f->num
        ){
            if (top->wall[fn]->exists
            ){
                PWall(p)->exists = TRUE;
              }
              else
              {
                PWall(p)->exists = FALSE;
            }
            if (top->edge[en]->exists
            ){
                PEdge(p)->exists = TRUE;
              }
              else
              {
                PEdge(p)->exists = FALSE;
            }
            PWall(p)->root = top->wall[fn]->root;
            PEdge(p)->root = top->edge[en]->root;
          }
        };


        /* Subdivision of degenerative @{edge->?}s */
        fprintf(stderr, "Eliminating degenerative @{edge->?}s\n");
        for (j = 0; j < top->NE; j++)
        {
          ??? e = top->edge[j];
          {
            if (e.degenerate)
            {
              ??? a = e.pa;
                ??? dfr = Octf.DegreeOfWall(a);
              {
                Subdivide@{Edge->?}(a,dfr,top);
              }
            }
          }
        }
      }
      
      if (o->element == Element.Wall)
      {
        /* Rescue the overall attributes of the triangulation */
        for (i = 0; i < top->wall.nelE; i++)
        {
          ??? p = top->wedge[i];
          Edge_t e = PEdge(p);
          with (en == e->num,
          double f = PWall(p), fn == f->num
        ){
            PWall(p)->root = top->wall[fn]->root;
            PEdge(p)->root = top->edge[en]->root;
          }
        };
    
        /* Subdivision of degenerative walls */
        fprintf(stderr, "Eliminating degenerative walls\n");

        for (j = 0; j < top->wall.nel; j++)
        {
          ??? f = top->wall[j];
            ??? a = f.pa;
          {
            if (f.degenerate)
            {
              SubdivideWall(a,top);
            }
          }
        }
      }

      if (o->element == Element.Tetrahedron)
      {
      /* Rescue the overall attributes of the triangulation */
        for (i = 0; i < top->wall.nelE; i++)
        {
          ??? p = top->wedge[i];
          Edge_t e = PEdge(p);
          with (en == e->num,
          double f = PWall(p), fn == f->num
        ){
            PWall(p)->root = top->wall[fn]->root;
            PEdge(p)->root = top->edge[en]->root;
          }
        };

        /* Subdivision of degenerative tetrahedra */
        fprintf(stderr, "Eliminating degenerative tetrahedra\n");
      
        for (j = 0; j < top->cell.nel; j++)
        {
          ??? p = top->cell[j];
            ??? r = Srot(top->cell[j]);
          Place_t a = Tors(r);
          {
            if (p.degenerate)
            {
              SubdivideTetrahedron(a,top);
            }
          }
        }
      }


      ??? newtop = MakeElemTable(top->wedge[0]);
      ??? nc = GenCoords(newtop)^;
      ??? cmt = "Created by Selective Subdivision of " & o->elename & " : " \
                   o->outFile;
      {
        WriteTopology(o->outFile, newtop, cmt & ".tp\n"
                                    & "on " & Today());

        WriteState(o->outFile,newtop,nc, cmt & ".st\n"
                                    & "on "&Today() & "\nRandom Geometry");
        
        WriteMaterials(o->outFile, newtop, cmt & ".ma\n"
                                    & "on " & Today());
      }      
    
    return 0;
  }


PROCEDURE Subdivide@{Edge->?}(
    Place_t @pn;
    uint n;
   ElemTableRec_t *top;
   ) ==
  VAR
    @{Node->?} x;
    a,bn,m1,m2,m3,t0,t1,Place_vec_t t2;
    wn,p: REF ARRAY OF @{Node->?};
  {

     a = NEW(REF Place_vec_t,n);
    bn = NEW(REF Place_vec_t,n); 
    m1 = NEW(REF Place_vec_t,n);
    m2 = NEW(REF Place_vec_t,n);
    m3 = NEW(REF Place_vec_t,n);
    t0 = NEW(REF Place_vec_t,n);
    t1 = NEW(REF Place_vec_t,n);
    t2 = NEW(REF Place_vec_t,n);

    /* save the @places */
    
    a = NEW(REF Place_vec_t,n);

    a[0] = an;
    ??? f = PWedge(a[0]).wall;
    ??? e = PEdge(a[0]);
    ??? ee = PEdge(NextE(a[0]));
    ??? ee_ = PEdge(PrevE(a[0]));
    {
      if (top->wall[f->num]->exists
      ){
          f->exists = TRUE;
        }
        else
        {
          f->exists = FALSE;
      }
      f.color  = top->wall[f->num].color;
      f.transp = top->wall[f->num].transp;
      f->root   = top->wall[f->num]->root;

      if (top->edge[e->num]->exists
      ){
          e->exists = TRUE;
        }
        else
        {
          e->exists = FALSE;
      }
      e.color  = top->edge[e->num].color;
      e.transp = top->edge[e->num].transp;
      e.radius = top->edge[e->num].radius;
      e->root   = top->edge[e->num]->root;

      if (top->edge[ee->num]->exists
      ){
           ee->exists = TRUE;
        }
        else
        {
           ee->exists = FALSE;
      }
      ee.color =  top->edge[ee->num].color;
      ee.transp = top->edge[ee->num].transp;
      ee.radius = top->edge[ee->num].radius;
      ee->root   = top->edge[e->num]->root;

      if (top->edge[ee_->num]->exists
      ){
          ee_->exists = TRUE;
        }
        else
        {
          ee_->exists = FALSE;
      }
      ee_.color =  top->edge[ee_->num].color;
      ee_.transp = top->edge[ee_->num].transp;
      ee_.radius = top->edge[ee_->num].radius;
      ee_->root   = top->edge[ee_->num]->root;
    }

    for (i = 1; i < n; i++)
    {
      a[i] = NextF(a[i-1]);
      ??? f = PWedge(a[i]).wall;
      ??? ee = PEdge(NextE(a[i]));
      ??? ee_ = PEdge(PrevE(a[i]));
      {
        if (top->wall[f->num]->exists
        ){
            f->exists = TRUE;
          }
          else
          {
            f->exists = FALSE;
        }
        f.color  = top->wall[f->num].color;
        f.transp = top->wall[f->num].transp;
        f->root   = top->wall[f->num]->root;

        if (top->edge[ee->num]->exists
         ){
             ee->exists = TRUE;
          }
          else
          {
             ee->exists = FALSE;
        }
        ee.color =  top->edge[ee->num].color;
        ee.transp = top->edge[ee->num].transp;
        ee.radius = top->edge[ee->num].radius;
        ee->root   = top->edge[ee->num]->root;

        if (top->edge[ee_->num]->exists
         ){
             ee_->exists = TRUE;
          }
          else
          {
             ee_->exists = FALSE;
        }
        ee_.color =  top->edge[ee_->num].color;
        ee_.transp = top->edge[ee_->num].transp;
        ee_.radius = top->edge[ee_->num].radius;
        ee_->root   = top->edge[ee_->num]->root;
      }
    }

    /* save the nodes */
    wn = NEW(REF ARRAY OF @{Node->?},n);
    for (i = 0; i < n; i++)
    {
      wn[i] = OrgV(PrevE(a[i]));
    }

    /* save the tetrahedra */
    p = NEW(REF ARRAY OF @{Node->?},n);
    for (i = 0; i < n; i++)
    {
      p[i] = PnegP(a[i]);
    }

    /* save other @places */
    bn = NEW(REF Place_vec_t,n);
    for (i = 0; i < n; i++)
    {
      bn[i] = Clock(PrevE(NextF(PrevE(a[i]))));
      ??? e =  PEdge(bn[i]);
        ??? f = PWedge( bn[i]).wall;
      {
        if (top->edge[e->num]->exists
        ){
            e->exists = TRUE;
          }
          else
          {
            e->exists = FALSE;
        }
        e.color  = top->edge[e->num].color;
        e.transp = top->edge[e->num].transp;
        e.radius = top->edge[e->num].radius;
        e->root   = top->edge[e->num]->root;

        if (top->wall[f->num]->exists
        ){
            f->exists = TRUE;
          }
          else
          {
            f->exists = FALSE;
        }
        f.color  = top->wall[f->num].color;
        f.transp = top->wall[f->num].transp;
        f->root   = top->wall[f->num]->root;
      }
    }

    /* insert new places */    
 
    for (i = 0; i < n; i++)
    {
      m1[i] = MakeWedge(); /* @{edge->?} component always FALSE */
      m2[i] = MakeWedge(); /* @{edge->?} component will depend  */
      m3[i] = MakeWedge(); /* @{edge->?} component always FALSE */
      t0[i] = MakeTriangle();  /* @{edge->?} component will depend  */
      t1[i] = NextE(t0[i]);    /* @{edge->?} component always FALSE */
      t2[i] = NextE(t1[i]);    /* @{edge->?} component always FALSE */
      ??? f0 = PWedge(t0[i]).wall;
        ??? f1 = PWedge(t1[i]).wall;
      ??? f2 = PWedge(t2[i]).wall;
      ??? e0 = PEdge(t0[i]);
      ??? e1 = PEdge(t1[i]);
      ??? e2 = PEdge(t2[i]);
      with (
       
      double em1 =  PEdge(m1[i]);
      double em2 =  PEdge(m2[i]);
      double em3 =  PEdge(m3[i]);
      double fm1 =  PWedge(m1[i]).wall;
      double fm2 =  PWedge(m2[i]).wall;
      double fm3 =  PWedge(m3[i]).wall
    ){
        /* with respect to m1[i], m2[i] and m3[i] */
        /* @{edge->?}s */
        em1->exists = FALSE;
        em3->exists = FALSE;
        if ((PEdge(a[i])->exists
        )){ em2->exists = TRUE;
          }
          else
          { em2->exists = FALSE;
        }
        em2->root =  PEdge(a[i])->root;      

        /* walls */
        if ((PWedge(a[i]).wall->exists
        )){
            fm1->exists = TRUE; 
            fm2->exists = TRUE;
            fm3->exists = TRUE;
          }
          else
          {
            fm1->exists = FALSE; 
            fm2->exists = FALSE;
            fm3->exists = FALSE;
        }

        /* with respect to insert wall */
        /* @{edge->?}s */
        e1->exists = FALSE;
        e2->exists = FALSE;
        if ((PEdge(bn[i])->exists
        )){ e0->exists = TRUE;
          }
          else
          { e0->exists = FALSE;
        }
        /* walls */
        f0->exists  = FALSE;
        f1->exists  = FALSE;
        f2->exists  = FALSE;
      }
    }

    /* Now, subdivide @{edge->?} and extend the subdivision on the edge's stars */
    
    x = MakeNode(0);
    ??? v = (Node_t)(x);
  {
  ){
      v.label = "VE";
      if ((PEdge(a[0])->exists))
      {
        v->exists = TRUE;
        v.radius = PEdge(a[0]).radius;
        v.color  = PEdge(a[0]).color;
        v.transp = PEdge(a[0]).transp;
      }
      else
      {
        v->exists = FALSE;
      }
      v->num    = NVE;
      NVE++;
    }
    for (j = 0; j < n; j++)
    {
      Place_t b = NextE(a[j]);
        Edge_t be = PEdge(b);
      Place_t c = NextE(b);
      Edge_t ce = PEdge(c);
      ??? u = OrgV(a[j]);
      ??? v = OrgV(b);
      ??? w = OrgV(c);
      with (
        /* save the attributes of the edge-wall component of the @place a[j] */
      double f = PWedge(a[j]).wall;
      double fn = f->num;
      double ff =  top->wall[fn];
      double fe = ff->exists;
      double fc = ff.color;
      double ft = ff.transp;
      double g = PWedge(m3[j]).wall;
      double ge = g->exists;
      double gc = g.color;
      double gt = g.transp;
      double p = PWedge(m1[j]).wall;
      double pe = p->exists;
      double pc = p.color;
      double pt = p.transp;
      double h = PEdge(m3[j])
    ){
        
        SetNextE(a[j],m1[j]);
        SetNextE(m1[j],c);
        SetNextE(m2[j],m3[j]);
        SetNextE(m3[j],Clock(b));
        
        SetOrg(a[j], u);
        SetOrg(Clock(a[j]), x);
        SetOrg(m2[j],v);
        SetOrg(Clock(m2[j]), x);
        SetOrg(m3[j], x);
        SetOrg(Clock(m3[j]), w);
        SetOrg(m1[j], x);
        SetOrg(Clock(m1[j]), w);
        SetNextF(m1[j],m3[j]);

        /* set the attributes for the wall component */
        
        SetRingWallInfo(a[j],f);
        SetRingWallInfo(m3[j],g);

        if (fe)
        {
          ge = TRUE;  pe = TRUE;
        }
        else
        {
          ge = FALSE; pe = FALSE;
        }
        gc = fc; pc = fc;
        gt = ft; pt = ft;

        SetRingEdgeInfo(m3[j],h);
        SetRingEdgeInfo(b,be);
        SetRingEdgeInfo(c,ce);
        SetRingWallInfo(bn[j],PWedge(bn[j]).wall);
        SetRingWallInfo(NextF(bn[j]),PWedge(NextF(bn[j])).wall);
      }
    }

    for (j = 0; j < n; j++)
    {
      SetNextF(Clock(m2[j]),Clock(m2[(j+1) % n]));
    }

    ??? k = PEdge(m2[0]);
    ??? q = PWedge(m2[0]).wall;
    ??? f = PWedge(a[0]).wall;
    ??? e = PEdge(a[0]);
    ??? en = e->num;
    ??? ea = top->edge[en];
    ??? ee = ea->exists;
    ??? er = ea.radius;
    ??? ec = ea.color;
    ??? et = e.transp;
    ??? fn = f->num;
    ??? fa = top->wall[fn];
    ??? fe = fa->exists;
    ??? fc = fa.color;
    ??? ft = fa.transp;
    ??? ke = k->exists;
    ??? kr = k.radius;
    ??? kc = k.color;
    ??? kt = k.transp;
    ??? qe = q->exists;
    ??? qc = q.color;
    ??? qt = q.transp;
    {
      if (ee){ ke = TRUE; }else{ ke = FALSE; }
      kr = er;
      kc = ec;
      kt = et;
      SetRingEdgeInfo(a[0],e);
      SetRingEdgeInfo(m2[0],k);
      if (fe){ qe = TRUE; }else{ qe = FALSE; }
      qc = fc;
      qt = ft;
    }

    for (j = 0; j < n; j++)
    {
      Place_t cn = NextF(bn[j]);
        ??? e1 = PEdge(t1[j]);
      ??? e2 = PEdge(t2[j]);
      Wall_t cnf = PWall(cn);
      Edge_t cne = PEdge(cn);
      {
        if (top->wall[cnf->num]->exists
        ){
            cnf->exists = TRUE;
          }
          else
          {
            cnf->exists = FALSE;
        }
        cnf.color  = top->wall[cnf->num].color;
        cnf.transp = top->wall[cnf->num].transp;

        if (top->edge[cne->num]->exists
        ){
            cne->exists = TRUE;
          }
          else
          {
            cne->exists = FALSE;
        }
        cne.color  = top->edge[cne->num].color;
        cne.transp = top->edge[cne->num].transp;
        cne.radius = top->edge[cne->num].radius;

        /* set the origins */

        SetOrgAll(t0[j],wn[j]);
        SetOrgAll(t1[j],wn[(j+1) % n]);
        SetOrgAll(t2[j],x);

        SetNextF(bn[j], t0[j]);
        SetNextF(t0[j],cn);

        SetNextF(m1[j],t2[j]);
        SetNextF(t2[j],m3[j]);

        SetNextF(Clock(m1[(j+1) % n]), t1[j]);
        SetNextF(t1[j],Clock(m3[(j+1) % n]));
        SetRingEdgeInfo(t0[j],PEdge(cn));
        SetRingEdgeInfo(t1[j],e1);
        SetRingEdgeInfo(t2[j],e2);
      }
    }    

    /* insert cells */
    for (j = 0; j < n; j++)
    {
      ??? q = Triangulation.MakeCell(0);
      {
        SetPnegOfNearbyWalls(a[j],p[j]);
        SetPnegOfNearbyWalls(m2[j],q);
      }
    }
  } /* END Subdivide@{Edge->?} */

PROCEDURE SubdivideWall(Place_t @p, ElemTableRec_t *top)
  /*
     Subdivide a degenerative triangular wall in three new walls through the
     insertion of a new medial node of type "VF" and six wall-@{edge->?} @places.
     Then, expand this subdivision on the wall's star. */
  VAR
    @{Node->?} x;
  {    
    /* Create the medial node "VF" */    
    x = MakeNode(0);
    ??? v = (Node_t)(x);
  {
  ){
      v->exists= FALSE;
      v.label = "VF";
      v->num   = NVF;
      NVF++;
    }      
    
    Wall_t fa = PWall(a);
    Edge_t ea = PEdge(a);
    Place_t b = NextE(a);
    Edge_t eb = PEdge(b);
    Place_t c = PrevE(a);
    Edge_t ec = PEdge(c);
    Place_t af = NextF(a);
    Place_t bf = NextF(b);
    Place_t cf = NextF(c);
    Place_t af_ = PrevF(a);
    Place_t bf_ = PrevF(b);
    Place_t cf_ = PrevF(c);
    Edge_t eaf = PEdge(PrevE(af));
    Edge_t eaf_ = PEdge(PrevE(af_));
    Wall_t faf = PWall(af);
    Wall_t faf_ = PWall(af_);
    Edge_t ebf = PEdge(PrevE(bf));
    Edge_t ebf_ = PEdge(PrevE(bf_));
    Wall_t fbf = PWall(bf);
    Wall_t fbf_ = PWall(bf_);
    Edge_t ecf = PEdge(PrevE(cf));
    Edge_t ecf_ = PEdge(PrevE(cf_));
    Wall_t fcf = PWall(cf);
    Wall_t fcf_ = PWall(cf_);
    ??? u = OrgV(a);
    ??? v = OrgV(b);
    ??? w = OrgV(c);
    Place_t d = MakeWedge();
    Place_t e = MakeWedge();
    Place_t f = MakeWedge();
    Place_t g = MakeWedge();
    Place_t h = MakeWedge();
    Place_t i = MakeWedge();
    Wall_t df = PWall(d);
    Wall_t ff = PWall(f);
    Wall_t hf = PWall(h);
    with (
      /* new walls to insert on the wall's star */
      f1  == MakeTriangle(),
      f2  == MakeTriangle(),
      f3  == MakeTriangle(),
      f4  == MakeTriangle(),
      f5  == MakeTriangle(),
      f6  == MakeTriangle(),
    double q1 = MakeCell(0);
    double q2 = MakeCell(1);
    double q3 = MakeCell(2);
    double q4 = MakeCell(3);
    double q5 = MakeCell(4);
    double q6 = MakeCell(5)
  ){
      ??? t = top->wall[fa->num];
      ??? tc = t.color;
      ??? tt = t.transp;
      ??? tr = t->root;
      {
        if (t->exists)
        {
          df->exists = TRUE;
          ff->exists = TRUE;
          hf->exists = TRUE;
          df.color  = tc;
          df.transp = tt;
          df->root   = tr;
          ff.color  = tc;
          ff.transp = tt;
          ff->root   = tr;
          hf.color  = tc;
          hf.transp = tt;
          hf->root   = tr;
        }
        else
        {
          df->exists = FALSE;
          ff->exists = FALSE;
          hf->exists = FALSE;
          df->root   = tr;
          ff->root   = tr;
          hf->root   = tr;
        }
      }

      /* first link connecting the @place "a" */

      SetNextE(a,d);
      SetNextE(d,e);
      SetNextE(e,a);
      SetRingWallInfo(d,df);

      /* second link connecting the @place "b" */
      
      SetNextE(b,f);
      SetNextE(f,g);
      SetNextE(g,b);
      SetRingWallInfo(f,ff);
      
      /* third link connecting the @place "c" */
      
      SetNextE(c,h);
      SetNextE(h,i);
      SetNextE(i,c);
      SetRingWallInfo(h,hf);

      /* save information about the original edges "PEdge(a)":  "ea" */

      ea->exists = top->edge[ea->num]->exists;
      ea.color  = top->edge[ea->num].color;
      ea.transp = top->edge[ea->num].transp;
      ea.radius = top->edge[ea->num].radius;

      
      /* save information about the original edges "PEdge(b)":  "eb" */

      eb->exists = top->edge[eb->num]->exists;
      eb.color  = top->edge[eb->num].color;
      eb.transp = top->edge[eb->num].transp;
      eb.radius = top->edge[eb->num].radius;
      

      /* save information about the original edges "PEdge(c)":  "ec" */

      ec->exists = top->edge[ec->num]->exists;
      ec.color  = top->edge[ec->num].color;
      ec.transp = top->edge[ec->num].transp;
      ec.radius = top->edge[ec->num].radius;

      /* save information about the wall "PWall(af)" : "faf" */

      faf->exists = top->wall[faf->num]->exists;
      faf.color  = top->wall[faf->num].color;
      faf.transp = top->wall[faf->num].transp;

      /* save information about the edge "PEdge(af)" : "eaf" */
      
      eaf->exists = top->edge[eaf->num]->exists;
      eaf.color  = top->edge[eaf->num].color;
      eaf.transp = top->edge[eaf->num].transp;
      eaf.radius = top->edge[eaf->num].radius;

      /* save information about the wall "PWall(af_)" : "faf_" */

      faf_->exists = top->wall[faf_->num]->exists;
      faf_.color = top->wall[faf_->num].color;
      faf_.transp= top->wall[faf_->num].transp;

      /* save information about the edge "PEdge(af_)" : "eaf_" */
      
      eaf_->exists = top->edge[eaf_->num]->exists;
      eaf_.color  = top->edge[eaf_->num].color;
      eaf_.transp = top->edge[eaf_->num].transp;
      eaf_.radius = top->edge[eaf_->num].radius;

      /* save information about the wall "PWall(bf)" : "fbf" */

      fbf->exists = top->wall[fbf->num]->exists;
      fbf.color = top->wall[fbf->num].color;
      fbf.transp= top->wall[fbf->num].transp;

      /* save information about the edge "PEdge(bf)" : "ebf" */

      ebf->exists = top->edge[ebf->num]->exists;
      ebf.color  = top->edge[ebf->num].color;
      ebf.transp = top->edge[ebf->num].transp;
      ebf.radius = top->edge[ebf->num].radius;

      /* save information about the wall "PWall(bf_)" : "fbf_" */

      fbf_->exists= top->wall[fbf_->num]->exists;
      fbf_.color = top->wall[fbf_->num].color;
      fbf_.transp= top->wall[fbf_->num].transp;

      /* save information about the edge "PEdge(bf_)" : "ebf_" */

      ebf_->exists = top->edge[ebf_->num]->exists;
      ebf_.color  = top->edge[ebf_->num].color;
      ebf_.transp = top->edge[ebf_->num].transp;
      ebf_.radius = top->edge[ebf_->num].radius;

      /* save information about the wall "PWall(cf)" : "fcf" */

      fcf->exists = top->wall[fcf->num]->exists;
      fcf.color = top->wall[fcf->num].color;
      fcf.transp= top->wall[fcf->num].transp;

      /* save information about the edge "PEdge(cf)" : "ecf" */

      ecf->exists = top->edge[ecf->num]->exists;
      ecf.color  = top->edge[ecf->num].color;
      ecf.transp = top->edge[ecf->num].transp;
      ecf.radius = top->edge[ecf->num].radius;

      /* save information about the wall "PWall(cf_)" : "fcf_" */

      fcf_->exists = top->wall[fcf_->num]->exists;
      fcf_.color = top->wall[fcf_->num].color;
      fcf_.transp= top->wall[fcf_->num].transp;

      /* save information about the edge "PEdge(cf_)" : "ecf_" */

      ecf_->exists = top->edge[ecf_->num]->exists;
      ecf_.color  = top->edge[ecf_->num].color;
      ecf_.transp = top->edge[ecf_->num].transp;
      ecf_.radius = top->edge[ecf_->num].radius;

      /* set the attributes for the internal walls */

      PWall(f1)->exists = FALSE;
      PWall(f2)->exists = FALSE;
      PWall(f3)->exists = FALSE;
      PWall(f4)->exists = FALSE;
      PWall(f5)->exists = FALSE;
      PWall(f6)->exists = FALSE;

      /* Now make the connections in the interior of wall */

      SetNextF(g,Clock(d));
      SetNextF(i,Clock(f));
      SetNextF(e,Clock(h));

      /* Now subdivide the "PposP(a)" tetrahedron ( the superior tetrahedron
         in my plan */

      /* first we, insert f1 */
      
      SetNextF(d, Clock(f1));
      SetNextF(Clock(f1), Clock(g));
      assert(NextF(Clock(g)) == d);
      SetNextF(NextE(f1),NextE(af));
      SetNextF(Clock(PrevE(bf)), NextE(f1));

      /* now, we insert f2 */
      
      SetNextF(f, Clock(f2));
      SetNextF(Clock(f2), Clock(i));
      assert(NextF(Clock(i)) == f);
      SetNextF(NextE(f2),NextE(bf));
      SetNextF(Clock(PrevE(cf)), NextE(f2));

      /* now, we insert f3 */
      
      SetNextF(h, Clock(f3));
      SetNextF(Clock(f3), Clock(e));
      assert(NextF(Clock(e)) == h);
      SetNextF(NextE(f3),NextE(cf));
      SetNextF(Clock(PrevE(af)), NextE(f3));
          
      /* now make connections between the internal walls inserted */

      SetNextF(PrevE(f1),PrevE(f3));
      SetNextF(PrevE(f3),PrevE(f2));
      assert(NextF(PrevE(f2)) == PrevE(f1));

      /* set the superior axial @{edge->?} */

      SetRingEdgeInfo(PrevE(f1), PEdge(PrevE(f1))); /* OK */
      PEdge(PrevE(f1))->exists = FALSE;
      
      /* Now subdivide the "PnegP(a)" tetrahedron ( the inferior tetrahedron
         in my plan */

      /* now, we insert f4 */

      SetNextF(Clock(g), Clock(f4));
      assert(NextF(Clock(f4)) == d);
      
      SetNextF(NextE(af_), NextE(f4));
      SetNextF(NextE(f4), Clock(PrevE(bf_)));
      
      /* now, we insert f5 */

      SetNextF(Clock(i), Clock(f5));
      assert(NextF(Clock(f5)) == f);
      SetNextF(NextE(bf_), NextE(f5));
      SetNextF(NextE(f5), Clock(PrevE(cf_)));
 
      /* now, we insert f6 */

      SetNextF(Clock(e), Clock(f6));
      assert(NextF(Clock(f6)) == h);      
      SetNextF(NextE(cf_), NextE(f6));
      SetNextF(NextE(f6), Clock(PrevE(af_)));
       
      /* now make connections between the internal walls inserted */

      SetNextF(PrevE(f4),PrevE(f5));
      SetNextF(PrevE(f5),PrevE(f6));
      assert(NextF(PrevE(f6)) == PrevE(f4));

      /* set the inferior axial @{edge->?} */
      SetRingEdgeInfo(PrevE(f4), PEdge(PrevE(f4))); /* OK */
      PEdge(PrevE(f4))->exists = FALSE;
      
      /* set the origins */
     
      SetOrgAll(f1, x);   
      SetOrgAll(a, u);
      SetOrgAll(b, v);
      SetOrgAll(c, w);
      SetOrgAll(PrevE(f1), OrgV(PrevE(af)));
      SetOrgAll(PrevE(f4), OrgV(PrevE(af_)));

      /* set the new internal @{edge->?}s */
      SetRingEdgeInfo(f1,PEdge(f1));
      PEdge(f1)->exists = FALSE;
      
      SetRingEdgeInfo(f2,PEdge(f2));
      PEdge(f2)->exists = FALSE;
      
      SetRingEdgeInfo(f3,PEdge(f3));
      PEdge(f3)->exists = FALSE;
      
      /* set the original edgess and walls in the superior level */
      SetRingEdgeInfo(PrevE(af),eaf);
      SetRingEdgeInfo(PrevE(bf),ebf);
      SetRingEdgeInfo(PrevE(cf),ecf);
      
      SetRingWallInfo(af,faf);
      SetRingWallInfo(bf,fbf);
      SetRingWallInfo(cf,fcf);

      /* set the original edgess and walls in the inferior label */
      SetRingEdgeInfo(PrevE(af_),eaf_);      
      SetRingEdgeInfo(PrevE(bf_),ebf_);
      SetRingEdgeInfo(PrevE(cf_),ecf_);
      
      SetRingWallInfo(af_,faf_);
      SetRingWallInfo(bf_,fbf_);
      SetRingWallInfo(cf_,fcf_);

      /* set the original edgess in the half label */

      SetRingEdgeInfo(a,ea);
      SetRingEdgeInfo(b,eb);
      SetRingEdgeInfo(c,ec);

      SetPnegOfNearbyWalls(Clock(PrevE(f1)),q1);
      SetPnegOfNearbyWalls(Clock(PrevE(f2)),q2);
      SetPnegOfNearbyWalls(Clock(PrevE(f3)),q3);
      SetPnegOfNearbyWalls(PrevE(f4),q4);
      SetPnegOfNearbyWalls(PrevE(f5),q5);
      SetPnegOfNearbyWalls(PrevE(f6),q6);
    }
  } /* END SubdivideWall */


PROCEDURE SubdivideTetrahedron(Place_t @p, ElemTableRec_t *top)
  /*
     Subdivide a degenerative tetrahedron in four new tetrahedra through the
     insertion of a new medial node of type "VP" more six walls and three
     @{edge->?}s. */
  @{Node->?} *y;
  { /* Create the medial node "VP" */    
    y = MakeNode(NVP);
    ??? v = (Node_t)(y);
  {
  ){
      v->exists= FALSE;
      v.label = "VP";
      NVP++;
    }

    Place_t b = PrevF(a);
    Place_t c = PrevE(b);
    Place_t d = NextF(c);
    Place_t e = PrevE(a);
    Place_t f = PrevF(e);
    Place_t g = PrevE(f);
    Place_t h = NextF(g);
    Wall_t fa = PWall(a);
    Wall_t fb = PWall(b);
    Wall_t fg = PWall(g);
    Wall_t fh = PWall(h);
    Edge_t ea = PEdge(a);
    Edge_t ec = PEdge(c);
    Edge_t ee = PEdge(e);
    Edge_t eh = PEdge(h);
    Edge_t eeb = PEdge(NextE(b));
    Edge_t eea = PEdge(NextE(a));
    Place_t i = NextE(a);
    Place_t j = PrevF(i);
    Place_t k = NextE(b);
    Place_t l = NextF(k);
    ??? u = OrgV(a);
    ??? v = OrgV(Clock(a));
    ??? w = OrgV(c);
    ??? x = OrgV(e);
    ??? f1 = MakeTriangle();
    ??? f2 = MakeTriangle();
    ??? f3 = MakeTriangle();
    ??? f4 = MakeTriangle();
    ??? f5 = MakeTriangle();
    ??? f6 = MakeTriangle();
    Cell_t q1 = MakeCell(0);
    Cell_t q2 = MakeCell(1);
    Cell_t q3 = MakeCell(2);
    Cell_t q4 = MakeCell(3);
    {
      /* save attributes for the original walls: "fa", "fb", "fg", "fh". */ 
      
      fa->exists = top->wall[fa->num]->exists;
      fa.color  = top->wall[fa->num].color;
      fa.transp = top->wall[fa->num].transp;

      fb->exists = top->wall[fb->num]->exists;
      fb.color  = top->wall[fb->num].color;
      fb.transp = top->wall[fb->num].transp;

      fg->exists = top->wall[fg->num]->exists;
      fg.color  = top->wall[fg->num].color;
      fg.transp = top->wall[fg->num].transp;

      fh->exists = top->wall[fh->num]->exists;
      fh.color  = top->wall[fh->num].color;
      fh.transp = top->wall[fh->num].transp;
            
      /* save attributes for the original edgess: "ea", "ec", "ee", "eh", "eeb",
         and "eea".
      */ 

      ea->exists = top->edge[ea->num]->exists;
      ea.color  = top->edge[ea->num].color;
      ea.transp = top->edge[ea->num].transp;
      ea.radius = top->edge[ea->num].radius;
      
      ec->exists = top->edge[ec->num]->exists;
      ec.color  = top->edge[ec->num].color;
      ec.transp = top->edge[ec->num].transp;
      ec.radius = top->edge[ec->num].radius;
      
      
      ee->exists = top->edge[ee->num]->exists;
      ee.color  = top->edge[ee->num].color;
      ee.transp = top->edge[ee->num].transp;
      ee.radius = top->edge[ee->num].radius;
     
      eh->exists = top->edge[eh->num]->exists;
      eh.color  = top->edge[eh->num].color;
      eh.transp = top->edge[eh->num].transp;
      eh.radius = top->edge[eh->num].radius;
      

      eea->exists = top->edge[eea->num]->exists;
      eea.color  = top->edge[eea->num].color;
      eea.transp = top->edge[eea->num].transp;
      eea.radius = top->edge[eea->num].radius;
            

      eeb->exists = top->edge[eeb->num]->exists;
      eeb.color  = top->edge[eeb->num].color;
      eeb.transp = top->edge[eeb->num].transp;
      eeb.radius = top->edge[eeb->num].radius;
    
      /* set the attributes for the new internal walls */
      
      PWall(f1)->exists = FALSE;
      PWall(f2)->exists = FALSE;
      PWall(f3)->exists = FALSE;
      PWall(f4)->exists = FALSE;
      PWall(f5)->exists = FALSE;
      PWall(f6)->exists = FALSE;
      
      /* insert f1 */
      SetNextF(b,f1);
      SetNextF(f1,a);

      /* insert f2 */
      SetNextF(c,f2);
      SetNextF(f2,d);

      /* insert f3 */
      SetNextF(f,f3);
      SetNextF(f3,e);

      /* set the relations among f1,f2 and f3 */
      SetNextF(Clock(NextE(f2)),PrevE(f1));
      SetNextF(PrevE(f1),Clock(NextE(f3)));
      SetNextF(Clock(NextE(f3)),Clock(NextE(f2)));
      SetRingEdgeInfo(PrevE(f1), PEdge(PrevE(f1)));
      PEdge(PrevE(f1))->exists = FALSE;
      
      /* insert f4 */
      SetNextF(j,f4);
      SetNextF(f4,i);

      /* insert f5 */
      SetNextF(k,f5);
      SetNextF(f5,l);

      /* insert f6 */
      SetNextF(g,f6);
      SetNextF(f6,h);
      
      /* set the internal relations along @{edge->?} "yv" */

      SetNextF(PrevE(f5),PrevE(f4));
      SetNextF(PrevE(f4),Clock(NextE(f1)));
      SetNextF(Clock(NextE(f1)), PrevE(f5));
      SetRingEdgeInfo(Clock(NextE(f1)),PEdge(Clock(NextE(f1))));
      PEdge(NextE(f1))->exists = FALSE;
      
      /* set the internal relations along @{edge->?} "wy" */
      
      SetNextF(NextE(f5),Clock(PrevE(f6)));
      SetNextF(Clock(PrevE(f6)),Clock(PrevE(f2)));
      SetNextF(Clock(PrevE(f2)),NextE(f5));

      SetRingEdgeInfo(NextE(f5),PEdge(NextE(f5)));
      PEdge(NextE(f5))->exists = FALSE;
      
      /* set the internal relations along @{edge->?} "xy" */
       
      SetNextF(NextE(f6), NextE(f4));
      SetNextF(NextE(f4), Clock(PrevE(f3)));
      SetNextF(Clock(PrevE(f3)), NextE(f6));
      SetRingEdgeInfo(NextE(f4), PEdge(NextE(f4)));
      PEdge(NextE(f4))->exists = FALSE;
      
      /* set the overall @{edge->?} component */
      SetRingEdgeInfo(a, ea);
      SetRingEdgeInfo(c, ec);
      SetRingEdgeInfo(e, ee);
      SetRingEdgeInfo(i, eea);
      SetRingEdgeInfo(k, eeb);
      SetRingEdgeInfo(h, eh);

      SetRingWallInfo(a, fa);
      SetRingWallInfo(b, fb);
      SetRingWallInfo(g, fg);
      SetRingWallInfo(h, fh);

      /* set the origins */
      SetOrgAll(a,u);
      SetOrgAll(Clock(a),v);
      SetOrgAll(c,w);
      SetOrgAll(e,x);
      SetOrgAll(PrevE(f1),y);

      /* set the cells */
      SetPnegOfNearbyWalls(a,q1);
      SetPnegOfNearbyWalls(Spin(b),q2);
      SetPnegOfNearbyWalls(Spin(g),q3);
      SetPnegOfNearbyWalls(Spin(f6),q4);
    }
  } /* END SubdivideTetrahedron */


Options_t *GetOptions(int argc, char **argv)
  { Options_t *o = (Options_t *)malloc(sizeof(Options_t));
    argparser_t *pp = argparser_new(stderr, argc, argv);
    argparser_set_help(pp, PROG_NAME " version " PROG_VERS ", usage:\n" PROG_HELP);
    argparser_set_info(pp, PROG_INFO);
    argparser_process_help_info_options(pp);
    
    argparser_get_keyword(pp, "-inFile");                               
        o->inFile = argparser_get_next(pp);
        argparser_get_keyword(pp, "-outFile");                               
        o->outFile = argparser_get_next(pp);
        argparser_get_keyword(pp, "-element");  
        o->elename = argparser_get_next(pp);
        if (0 == strcmp(o->elename, "@{edge->?}")))
        {
          o->element = Element.edge
        }
        else
         if (0 == strcmp(o->elename, "wall"))){
          o->element = Element.Wall 
        }
        else
         if (0 == strcmp(o->elename, "tetrahedron"))){
          o->element = Element.Tetrahedron
        }
        else
        {
          argparser_error(pp, "bad element \"" & argparser_get_next(pp) & "\"\n")
        }
        argparser_finish(pp);                                       
      
  -----------------------------------
#define _HELP \
  fprintf(stderr, "Usage: SelectSubdivision\\\n" \
   "     -inFile <name> -outFile <name>\\\n" \
   "     -element { @{edge->?} | wall | tetrahedron }\n");
          
END¦                                              
      }
    }
    return o;
  } /* END GetOptions */

/* UNUSED */
PROCEDURE OldSubdivide@{Edge->?}(
  Place_t @pn;
  uint n;
 ElemTableRec_t *top;
 ) ==
        
    @{Node->?} x;
    a,bn,m1,m2,m3,t0,t1,Place_vec_t t2;
    wn,p: REF ARRAY OF @{Node->?};
  {  a = NEW(REF Place_vec_t,n);
    bn = NEW(REF Place_vec_t,n); 
    m1 = NEW(REF Place_vec_t,n);
    m2 = NEW(REF Place_vec_t,n);
    m3 = NEW(REF Place_vec_t,n);
    t0 = NEW(REF Place_vec_t,n);
    t1 = NEW(REF Place_vec_t,n);
    t2 = NEW(REF Place_vec_t,n);

    /* save the @places */
    
    a = NEW(REF Place_vec_t,n);
    a[0] = an;
    ??? f = PWedge(a[0]).wall;
    ??? e = PEdge(a[0]);
    ??? ee = PEdge(NextE(a[0]));
    ??? ee_ = PEdge(PrevE(a[0]));
    {
      if (top->wall[f->num]->exists
      ){
          f->exists = TRUE;
        }
        else
        {
          f->exists = FALSE;
      }
      f.color  = top->wall[f->num].color;
      f.transp = top->wall[f->num].transp;

      if (top->edge[e->num]->exists
      ){
          e->exists = TRUE;
        }
        else
        {
          e->exists = FALSE;
      }
      e.color  = top->edge[e->num].color;
      e.transp = top->edge[e->num].transp;
      e.radius = top->edge[e->num].radius;
      
      if (top->edge[ee->num]->exists
      ){
           ee->exists = TRUE;
        }
        else
        {
           ee->exists = FALSE;
      }
      ee.color =  top->edge[ee->num].color;
      ee.transp = top->edge[ee->num].transp;
      ee.radius = top->edge[ee->num].radius;

      if (top->edge[ee_->num]->exists
      ){
          ee_->exists = TRUE;
        }
        else
        {
          ee_->exists = FALSE;
      }
      ee_.color =  top->edge[ee_->num].color;
      ee_.transp = top->edge[ee_->num].transp;
      ee_.radius = top->edge[ee_->num].radius;
    }

    for (i = 1; i < n; i++)
    {
      a[i] = NextF(a[i-1]);
      ??? f = PWedge(a[i]).wall;
      ??? ee = PEdge(NextE(a[i]));
      ??? ee_ = PEdge(PrevE(a[i]));
      {
        if (top->wall[f->num]->exists
        ){
            f->exists = TRUE;
          }
          else
          {
            f->exists = FALSE;
        }
        f.color  = top->wall[f->num].color;
        f.transp = top->wall[f->num].transp;

        if (top->edge[ee->num]->exists
         ){
             ee->exists = TRUE;
          }
          else
          {
             ee->exists = FALSE;
        }
        ee.color =  top->edge[ee->num].color;
        ee.transp = top->edge[ee->num].transp;
        ee.radius = top->edge[ee->num].radius;

        if (top->edge[ee_->num]->exists
         ){
             ee_->exists = TRUE;
          }
          else
          {
             ee_->exists = FALSE;
        }
        ee_.color =  top->edge[ee_->num].color;
        ee_.transp = top->edge[ee_->num].transp;
        ee_.radius = top->edge[ee_->num].radius;
      }
    }

    /* save the nodes */
    wn = NEW(REF ARRAY OF @{Node->?},n);
    for (i = 0; i < n; i++)
    {
      wn[i] = OrgV(PrevE(a[i]));
    }

    /* save the tetrahedra */
    p = NEW(REF ARRAY OF @{Node->?},n);
    for (i = 0; i < n; i++)
    {
      p[i] = PnegP(a[i]);
    }

    /* save other @places */
    bn = NEW(REF Place_vec_t,n);
    for (i = 0; i < n; i++)
    {
      bn[i] = Clock(PrevE(NextF(PrevE(a[i]))));
      ??? e =  PEdge(bn[i]);
        ??? f =  PWedge(bn[i]).wall;
      {
        if (top->edge[e->num]->exists
        ){
            e->exists = TRUE;
          }
          else
          {
            e->exists = FALSE;
        }
        e.color  = top->edge[e->num].color;
        e.transp = top->edge[e->num].transp;
        e.radius = top->edge[e->num].radius;
        
        if (top->wall[f->num]->exists
        ){
            f->exists = TRUE;
          }
          else
          {
            f->exists = FALSE;
        }
        f.color  = top->wall[f->num].color;
        f.transp = top->wall[f->num].transp;
      }
    }

    /* insert wedges and edges */    
 
    for (i = 0; i < n; i++)
    {
      m1[i] = MakeWedge();
      m2[i] = MakeWedge();
      m3[i] = MakeWedge();
      t0[i] = MakeTriangle();
      t1[i] = NextE(t0[i]);
      t2[i] = NextE(t1[i]);
      ??? f = PWedge(t0[i]).wall;
        ??? e0 = PEdge(t0[i]);
      ??? e1 = PEdge(t1[i]);
      ??? e2 = PEdge(t2[i]);
      with (
        /* new tests */
      double em1i =  PEdge(m1[i]);
      double em2i =  PEdge(m2[i])
    ){
        /* new tests */
        
        em1i->exists = FALSE;
        em2i->exists = FALSE;
       
        /* end new tests */
        f->exists  = FALSE;
        f.transp  = (frgb_t){1.0,1.0,1.0};
        e1->exists = FALSE;
        e2->exists = FALSE;
        e1.radius = 0.003;
        e1.transp = (frgb_t){1.0,1.0,1.0};
        e2.radius = 0.003;
        e2.transp = (frgb_t){1.0,1.0,1.0};
        if ((PEdge(bn[i])->exists
        )){ e0->exists = TRUE;
          }
          else
          { e0->exists = FALSE;
        }
      }
    }

    /* Now, subdivide @{edge->?} and extend the subdivision on the edge's stars */
    
    x = MakeNode(0);
    ??? v = (Node_t)(x);
  {
  ){
      v.label = "VE";
      if ((PEdge(a[0])->exists))
      {
        v.radius = PEdge(a[0]).radius;
        v.color  = PEdge(a[0]).color;
        v.transp = PEdge(a[0]).transp;
      }
      else
      {
        v.radius = 0.00;
        v.color  = (frgb_t){1.0,1.0,1.0};
        v.transp = (frgb_t){1.0,1.0,1.0};
      }
      v->num    = NVE;
      NVE++;
    }
    for (j = 0; j < n; j++)
    {
      Place_t b = NextE(a[j]);
        Edge_t be = PEdge(b);
      Place_t c = NextE(b);
      Edge_t ce = PEdge(c);
      ??? u = OrgV(a[j]);
      ??? v = OrgV(b);
      ??? w = OrgV(c);
      with (
        /* save the attributes of the edge-wall component of the @place a[j] */
      double f = PWedge(a[j]).wall;
      double fe = f->exists;
      double fc = f.color;
      double ft = f.transp;
      double g = PWedge(m3[j]).wall;
      double ge = g->exists;
      double gc = g.color;
      double gt = g.transp;
      double h = PEdge(m3[j])
    ){
        SetNextE(a[j],m1[j]);
        SetNextE(m1[j],c);
        SetNextE(m2[j],m3[j]);
        SetNextE(m3[j],Clock(b));
        
        SetOrg(a[j], u);
        SetOrg(Clock(a[j]), x);
        SetOrg(m2[j],v);
        SetOrg(Clock(m2[j]), x);
        SetOrg(m3[j], x);
        SetOrg(Clock(m3[j]), w);
        SetOrg(m1[j], x);
        SetOrg(Clock(m1[j]), w);
        SetNextF(m1[j],m3[j]);

        /* set the attributes for the wall component */
        
        SetRingWallInfo(a[j],f);
        SetRingWallInfo(m3[j],g);

        if (fe)
        {
          ge = TRUE;
        }
        else
        {
          ge = FALSE;
        }
        gc = fc;
        gt = ft;

        SetRingEdgeInfo(m3[j],h);
        SetRingEdgeInfo(b,be);
        SetRingEdgeInfo(c,ce);
        SetRingWallInfo(PWedge(bn[j],bn[j]).wall);
        SetRingWallInfo(NextF(bn[j]),PWedge(NextF(bn[j])).wall);
      }
    }

    for (j = 0; j < n; j++)
    {
      SetNextF(Clock(m2[j]),Clock(m2[(j+1) % n]));
    }

    ??? k = PEdge(m2[0]);
    ??? e = PEdge(a[0]);
    ??? ee = e->exists;
    ??? er = e.radius;
    ??? ec = e.color;
    ??? et = e.transp;
    ??? ke = k->exists;
    ??? kr = k.radius;
    ??? kc = k.color;
    ??? kt = k.transp;
    {
      if (ee){ ke = TRUE; }else{ ke = FALSE; }
      kr = er;
      kc = ec;
      kt = et;
      SetRingEdgeInfo(a[0],e);
      SetRingEdgeInfo(m2[0],k);
    }

    for (j = 0; j < n; j++)
    {
      Place_t cn = NextF(bn[j]);
        ??? e1 = PEdge(t1[j]);
      ??? e2 = PEdge(t2[j]);
      Wall_t cnf = PWall(cn);
      Edge_t cne = PEdge(cn);
      {
        if (top->wall[cnf->num]->exists
        ){
            cnf->exists = TRUE;
          }
          else
          {
            cnf->exists = FALSE;
        }
        cnf.color  = top->wall[cnf->num].color;
        cnf.transp = top->wall[cnf->num].transp;

        if (top->edge[cne->num]->exists
        ){
            cne->exists = TRUE;
          }
          else
          {
            cne->exists = FALSE;
        }
        cne.color  = top->edge[cne->num].color;
        cne.transp = top->edge[cne->num].transp;
        cne.radius = top->edge[cne->num].radius;
        
        /* set the origins */

        SetOrgAll(t0[j],wn[j]);
        SetOrgAll(t1[j],wn[(j+1) % n]);
        SetOrgAll(t2[j],x);

        SetNextF(bn[j], t0[j]);
        SetNextF(t0[j],cn);

        SetNextF(m1[j],t2[j]);
        SetNextF(t2[j],m3[j]);

        SetNextF(Clock(m1[(j+1) % n]), t1[j]);
        SetNextF(t1[j],Clock(m3[(j+1) % n]));
        SetRingEdgeInfo(t0[j],PEdge(cn));
        SetRingEdgeInfo(t1[j],e1);
        SetRingEdgeInfo(t2[j],e2);
      }
    }    

    /* insert cells */
    for (j = 0; j < n; j++)
    {
      ??? q = Triangulation.MakeCell(0);
      {
        SetPnegOfNearbyWalls(a[j],p[j]);
        SetPnegOfNearbyWalls(m2[j],q);
      }
    }
  } /* END OldSubdivide@{Edge->?} */

/* UNUSED */
uint Incident(
   Place_t b;
   uint n;
  ElemTableRec_t *top;
  )
  /* This module indicates the number of existing walls incident to
     the component @{edge->?} of the @place b. */
  VAR
     Place_t bn = b;
    in: uint = 0;
  { for (j = 0; j < n; j++){
      ??? fn = PWall(bn)->num;
        ??? f = top->wall[fn];
      {
        if (f->exists){ in++; }
        bn = NextF(bn);
      }
    }
    assert(bn == b);
    return in;
  } /* END Incident */


{
  DoIt()
} SelectSubdivision.


  /* Copyright © 2000     Universidade Estadual de Campinas (UNICAMP)     */