geompack_insfac.cxx

hl2 source code. Do not use it illegal.

/* insfac.f -- translated by f2c (version 19990311).
*/


#include <ivp_physics.hxx>
#include <ivp_betterdebugmanager.hxx>
#include <geompack.hxx>


void IVP_Geompack::insfac_(
			   int		*p,
			   double	*nrmlc,
			   int		*nce,
			   int		*cedge,
			   double	*cdang
			  ) {

    double *vcl                       = this->g_vcl;
    int	   *polyhedronfirstfaceoffset = this->g_polyhedronfirstfaceoffset;
    int	   *polyhedronfaceindices     = this->g_polyhedronfaceindices;
    int    *facesdata                 = this->g_facesdata;
    int    *facestype                 = this->g_facestype;
    double *normals                   = this->g_normals;
    int    *faceverticeslist          = this->g_faceverticeslist;
    double *edge_angles               = this->g_edge_angles;

    /* System generated locals */
    int i__1, i__2;
    double d__1;

    /* Local variables */
    int head;
    long int docf;
    int pind, a, b, c__, f, g, i__, j, k, tailn, tailp, tailt;
    int la, lb;
    int fp, sf, nv, sp, np1;
    double ang;
    long int dof;
    int ptr;


/*     Written and copyright by: */
/*        Barry Joe, Dept. of Computing Science, Univ. of Alberta */
/*        Edmonton, Alberta, Canada  T6G 2H1 */
/*        Phone: [phone]      Email: [e-mail] */

/*     Purpose: Insert a new face (cut face) in polyhedral decomposition */
/*        data structure. It is assumed that interior of face does not */
/*        intersect any other faces. */

/*     Input parameters: */
/*        P - polyhedron index */
/*        NRMLC(1:3) - unit normal vector of cut plane */
/*        NCE - number of edges in cut face */
/*        CEDGE(1:2,0:NCE) - CEDGE(1,I) is an index of VCL, indices > ADDR_OF_N_ORIGINAL_VERTICES */
/*              are new points; CEDGE(2,I) = J indicates that edge of cut */
/*              face ending at CEDGE(1,I) is edge from J to FACEVERTICESLIST(SUCC,J) */
/*              if J > 0; else if J < 0 then edge of cut face ending at */
/*              CEDGE(1,I) is a new edge and CEDGE(1,I) lies on edge from */
/*              -J to FACEVERTICESLIST(SUC,-J) and new edge lies in face FACEVERTICESLIST(FACN,-J); */
/*              CEDGE(2,I) always refers to an edge in the subpolyhedron */
/*              in negative half-space; CEDGE(1,NCE) = CEDGE(1,0); */
/*              CEDGE(2,0) is not input but is used temporarily */
/*        CDANG(1:NCE) - dihedral angles created by edges of cut polygon */
/*              in pos. half-space; neg. sign for angle I indicates that */
/*              face containing edge I is oriented CW in polyhedron P */
/*        ADDR_OF_N_ORIGINAL_VERTICES - number of vertex coordinates (excluding new ones)) */
/*        NFACE - number of faces or positions used in FACESDATA array */
/*        ADDR_OF_N_WORK_VERTICES - number of positions used in FACEVERTICESLIST, EDGE_ANGLES arrays */
/*        NPOLH - number of polyhedra or positions used in POLYHEDRONFIRSTFACEOFFSET array */
/*        ADDR_OF_N_POLYHEDRONFACES - number of positions used in POLYHEDRONFACEINDICES array */
/*        ADDR_OF_SIZE_FACEARRAYS - maximum size available for FACESDATA, FACESTYPE, NORMALS arrays */
/*        ADDR_OF_SIZE_FACEVERTEXARRAYS - maximum size available for FACEVERTICESLIST, EDGE_ANGLES arrays */
/*        ADDR_OF_SIZE_POLYHEDRONFIRSTFACEOFFSET - maximum size available for POLYHEDRONFIRSTFACEOFFSET array */
/*        ADDR_OF_SIZE_POLYHEDRONFACEINDICES - maximum size available for POLYHEDRONFACEINDICES array */
/*        VCL(1:3,1:ADDR_OF_N_ORIGINAL_VERTICES+?) - vertex coordinate list; the new vertices to */
/*           be inserted as indicated by CEDGE are after column ADDR_OF_N_ORIGINAL_VERTICES */
/*        FACESDATA(1:3,1:NFACE) - face pointer list */
/*        FACESTYPE(1:NFACE) - face types */
/*        NORMALS(1:3,1:NFACE) - unit normal vectors for faces */
/*        FACEVERTICESLIST(1:6,1:ADDR_OF_N_WORK_VERTICES) - face vertex list */
/*        EDGE_ANGLES(1:ADDR_OF_N_WORK_VERTICES) - edge angles */
/*        POLYHEDRONFIRSTFACEOFFSET(1:NPOLH) - head ptr to face indices in POLYHEDRONFACEINDICES for each polyh */
/*        POLYHEDRONFACEINDICES(1:2,1:ADDR_OF_N_POLYHEDRONFACES) - list of signed face indices for each polyh */

/*     Updated parameters: */
/*        CEDGE(1:1,1:NCE) - updated to edges of cut face wrt positive */
/*              half-space */
/*        CEDGE(2:2,1:NCE) - neg. entries updated to index of new edge */
/*        ADDR_OF_N_ORIGINAL_VERTICES,NFACE,ADDR_OF_N_WORK_VERTICES,NPOLH,ADDR_OF_N_POLYHEDRONFACES - sizes updated due to cut face */
/*        FACESDATA,FACESTYPE,NORMALS,FACEVERTICESLIST,EDGE_ANGLES,POLYHEDRONFIRSTFACEOFFSET,POLYHEDRONFACEINDICES - updated by cut face */

/*     Routines called: */
/*        INSVR3, INSED3 */

/*     Abnormal return: */
/*        IERR is set to 15, 16, 17, or 18 */




/*     Insert new vertices and update CEDGE(2,*). */

    /* Parameter adjustments */
    --nrmlc;
    --cdang;
    normals -= 4;
    --facestype;
    facesdata -= 4;
    --edge_angles;
    faceverticeslist -= 7;
    --polyhedronfirstfaceoffset;
    polyhedronfaceindices -= 3;

    /* Function Body */
    i__1 = *nce - 1;
    for (i__ = 0; i__ <= i__1; ++i__) {
	if (cedge[(i__ << 1) + 1] <= this->n_original_vertices ) {
	    goto L10;
	}
	if (i__ == 0) {
	    j = *nce;
	} else {
	    j = i__;
	}
	a = -cedge[(j << 1) + 2];

	// *******************************************************************
	insvr3_(a);

	faceverticeslist = this->g_faceverticeslist;
	faceverticeslist -= 7;

	edge_angles = this->g_edge_angles;
	--edge_angles;


	if (this->ierr != 0) {
	    return;
	}
	if (cdang[j] < 0.) {
	    cedge[(j << 1) + 2] = -faceverticeslist[a * 6 + 3];
	}
L10:
	;
    }

/*     Insert new edges and update CEDGE(2,*). */

    cedge[2] = cedge[(*nce << 1) + 2];
    i__1 = *nce;
    for (i__ = 1; i__ <= i__1; ++i__) {
	b = -cedge[(i__ << 1) + 2];
	if (b < 0) {
	    goto L30;
	}
	f = faceverticeslist[b * 6 + 2];
	la = cedge[((i__-1) << 1) + 1];
	a = cedge[((i__-1) << 1) + 2];
	if (a < 0) {
/*           This can only occur if I = 1. */
	    a = -a;
	    if (faceverticeslist[a * 6 + 1] == la) {
		a = faceverticeslist[a * 6 + 4];
	    } else {
		a = faceverticeslist[a * 6 + 3];
	    }
	}
L20:
	if (faceverticeslist[a * 6 + 1] == la) {
	    a = faceverticeslist[a * 6 + 4];
	    j = la - faceverticeslist[a * 6 + 1];
	    sf = *p;
	} else {
	    a = faceverticeslist[a * 6 + 3];
	    j = faceverticeslist[faceverticeslist[a * 6 + 3] * 6 + 1] - la;
	    sf = -(*p);
	}
	if (j * sf > 0) {
	    a = faceverticeslist[a * 6 + 6];
	} else {
	    a = faceverticeslist[a * 6 + 5];
	}
	fp = faceverticeslist[a * 6 + 2];
	if (fp != f) {
	    goto L20;
	}
	if (faceverticeslist[a * 6 + 1] == la) {
	    j = a;
	    a = faceverticeslist[b * 6 + 3];
	    b = j;
	} else {
	    a = faceverticeslist[a * 6 + 3];
	}

	// *******************************************************************
	insed3_(&a,
		&b,
		&polyhedronfirstfaceoffset[1]
		);

	polyhedronfaceindices = this->g_polyhedronfaceindices;
	polyhedronfaceindices -= 3;

	facesdata = this->g_facesdata;
	facesdata -= 4;

	facestype = this->g_facestype;
	--facestype;

	normals = this->g_normals;
	normals -= 4;

	faceverticeslist = this->g_faceverticeslist;
	faceverticeslist -= 7;

	edge_angles = this->g_edge_angles;
	--edge_angles;

	if (this->ierr != 0) {
	    return;
	}
	cedge[(i__ << 1) + 2] = a;
L30:
	;
    }

/*     Insert cut face into decomposition data structure. Subpolyhedron */
/*     in negative half space is numbered P, other is numbered NPOLH. */

    this->nface++;
    this->npolh++;
    this->n_polyhedronfaces += 2;
recheck_size:
    if ( (this->n_work_vertices + *nce + 2) > this->size_facevertexarrays) {
	int res = 0;
	res = res | increase_memory((void **)&this->g_faceverticeslist, &this->size_facevertexarrays, 6*sizeof(int));
	res = res | increase_memory((void **)&this->g_edge_angles     , &this->size_facevertexarrays,   sizeof(double));
	if ( res == 0 ) {
	    this->ierr = 500;
	    return;
	}
//	*addr_of_size_facevertexarrays *= 2;
	this->size_facevertexarrays +=1024;
	faceverticeslist = this->g_faceverticeslist;
	faceverticeslist -= 7;
	edge_angles = this->g_edge_angles;
	--edge_angles;
	goto recheck_size;
    }

    else if ( (this->nface+2) > this->size_facearrays ) {
	int res = 0;
	res = res | increase_memory((void **)&this->g_normals  , &this->size_facearrays, 3*sizeof(double));
	res = res | increase_memory((void **)&this->g_facesdata, &this->size_facearrays, 3*sizeof(int));
	res = res | increase_memory((void **)&this->g_facestype, &this->size_facearrays,   sizeof(int));
	if ( res == 0 ) {
	    this->ierr = 500;
	    return;
	}
//	*addr_of_size_facearrays *= 2;
	this->size_facearrays += 1024;
	facesdata = this->g_facesdata;
	facesdata -= 4;
	facestype = this->g_facestype;
	--facestype;
	normals = this->g_normals;
	normals -= 4;
	goto recheck_size;
    }

    else if ( (2+this->n_polyhedronfaces) > this->size_polyhedronfaceindices ) {
	int res = increase_memory((void **)&this->g_polyhedronfaceindices, &this->size_polyhedronfaceindices, 2*sizeof(int));
	if ( res == 0 ) {
	    this->ierr = 500;
	    return;
	}
//	*addr_of_size_polyhedronfaceindices *= 2;
	this->size_polyhedronfaceindices += 1024;
	polyhedronfaceindices = this->g_polyhedronfaceindices;
	polyhedronfaceindices -= 3;
	goto recheck_size;
    }

    else if ( (1+this->npolh) > this->size_polyhedronfirstfaceoffset ) {
	int res = increase_memory((void **)&this->g_polyhedronfirstfaceoffset, &this->size_polyhedronfirstfaceoffset, sizeof(int));
	if ( res == 0 ) {
	    this->ierr = 500;
	    return;
	}
//	*addr_of_size_polyhedronfirstfaceoffset *= 2;
	this->size_polyhedronfirstfaceoffset += 1024;
	polyhedronfirstfaceoffset = this->g_polyhedronfirstfaceoffset;
	--polyhedronfirstfaceoffset;
	goto recheck_size;
    }

    nv = this->n_work_vertices;
    facesdata[this->nface * 3 + 1] = this->n_work_vertices + 1;
    facesdata[this->nface * 3 + 2] = *p;
    facesdata[this->nface * 3 + 3] = -this->npolh;
    facestype[this->nface] = 0;
    normals[this->nface * 3 + 1] = nrmlc[1];
    normals[this->nface * 3 + 2] = nrmlc[2];
    normals[this->nface * 3 + 3] = nrmlc[3];
    i__1 = *nce - 1;
    for (i__ = 0; i__ <= i__1; ++i__) {
	this->n_work_vertices++;
	faceverticeslist[this->n_work_vertices * 6 + 1] = cedge[(i__ << 1) + 1];
	faceverticeslist[this->n_work_vertices * 6 + 2] = this->nface;
	faceverticeslist[this->n_work_vertices * 6 + 3] = this->n_work_vertices + 1;
	faceverticeslist[this->n_work_vertices * 6 + 4] = this->n_work_vertices - 1;
/* L40: */
    }
    faceverticeslist[this->n_work_vertices * 6 + 3] = facesdata[this->nface * 3 + 1];
    faceverticeslist[facesdata[this->nface * 3 + 1] * 6 + 4] = this->n_work_vertices;

/*     Set CEDGE(1,*) to edges of face in polyh NPOLH (after split). New */
/*     face is CCW from outside new P which contains edges of CEDGE(2,*). */

    i__1 = *nce;