void
nmg_visit_shell(struct shell *s, const struct nmg_visit_handlers *htab, void *state)


/* Handler's private state */
{
    struct faceuse *fu;
    struct loopuse *lu;
    struct edgeuse *eu;

    NMG_CK_SHELL(s);

    if (htab->bef_shell) htab->bef_shell((uint32_t *)s, state, 0);

    for (BU_LIST_FOR(fu, faceuse, &s->fu_hd)) {
	nmg_visit_faceuse(fu, htab, state);
    }
    for (BU_LIST_FOR(lu, loopuse, &s->lu_hd)) {
	nmg_visit_loopuse(lu, htab, state);
    }
    for (BU_LIST_FOR(eu, edgeuse, &s->eu_hd)) {
	nmg_visit_edgeuse(eu, htab, state);
    }
    if (s->vu_p) nmg_visit_vertexuse(s->vu_p, htab, state);
    if (htab->vis_shell_a && s->sa_p)
	htab->vis_shell_a((uint32_t *)s->sa_p, state, 0);

    if (htab->aft_shell) htab->aft_shell((uint32_t *)s, state, 1);
}

void
rt_vlist_export(struct bu_vls *vls, struct bu_list *hp, const char *name)
{
    register struct bn_vlist *vp;
    size_t nelem;
    size_t namelen;
    size_t nbytes;
    unsigned char *buf;
    unsigned char *bp;

    BU_CK_VLS(vls);

    /* Count number of element in the vlist */
    nelem = 0;
    for (BU_LIST_FOR(vp, bn_vlist, hp)) {
	nelem += vp->nused;
    }

    /* Build output buffer for binary transmission
     * nelem[4], String[n+1], cmds[nelem*1], pts[3*nelem*8]
     */
    namelen = strlen(name)+1;
    nbytes = namelen + 4 + nelem * (1+ELEMENTS_PER_VECT*SIZEOF_NETWORK_DOUBLE) + 2;

    /* FIXME: this is pretty much an abuse of vls.  should be using
     * vlb for variable-length byte buffers.
     */
    bu_vls_setlen(vls, (int)nbytes);
    buf = (unsigned char *)bu_vls_addr(vls);
    *(uint32_t *)buf = htonl((uint32_t)nelem);
    bp = buf+sizeof(uint32_t);
    bu_strlcpy((char *)bp, name, namelen);
    bp += namelen;

    /* Output cmds, as bytes */
    for (BU_LIST_FOR(vp, bn_vlist, hp)) {
	register int i;
	register int nused = vp->nused;
	register int *cmd = vp->cmd;
	for (i = 0; i < nused; i++) {
	    *bp++ = *cmd++;
	}
    }

    /* Output points, as three 8-byte doubles */
    for (BU_LIST_FOR(vp, bn_vlist, hp)) {
	register int i;
	register int nused = vp->nused;
	register point_t *pt = vp->pt;

	/* must be double for import and export */
	double point[ELEMENTS_PER_POINT];

	for (i = 0; i < nused; i++) {
	    VMOVE(point, pt[i]); /* convert fastf_t to double */
	    htond(bp, (unsigned char *)point, ELEMENTS_PER_VECT);
	    bp += ELEMENTS_PER_VECT*SIZEOF_NETWORK_DOUBLE;
	}
    }
}

extern "C" void
rt_nmg_brep(ON_Brep **b, const struct rt_db_internal *ip, const struct bn_tol *tol)
{
    struct model *m;
    struct nmgregion *r;
    struct shell *s;
    struct faceuse *fu;

    // Verify NMG
    RT_CK_DB_INTERNAL(ip);
    m = (struct model *)ip->idb_ptr;
    NMG_CK_MODEL(m);

    // Both NMG and brep structures re-use components between faces.  In order to track
    // when the conversion routine has already handled an NMG element, use an array.
    long *brepi = static_cast<long*>(bu_malloc(m->maxindex * sizeof(long), "rt_nmg_brep: brepi[]"));
    for (int i = 0; i < m->maxindex; i++) brepi[i] = -INT_MAX;

    // Iterate over all faces in the NMG
    for (BU_LIST_FOR(r, nmgregion, &m->r_hd)) {
	for (BU_LIST_FOR(s, shell, &r->s_hd)) {
	    for (BU_LIST_FOR(fu, faceuse, &s->fu_hd)) {
		NMG_CK_FACEUSE(fu);
		if (fu->orientation != OT_SAME) continue;
		if (nmg_brep_face(b, fu, tol, brepi)) return;
	    }
	    (*b)->SetTrimIsoFlags();
	}
    }
}

inline void
rt_metaball_norm_internal(vect_t *n, point_t *p, struct rt_metaball_internal *mb)
{
    struct wdb_metaballpt *mbpt;
    vect_t v;
    fastf_t a;

    VSETALL(*n, 0.0);

    switch (mb->method) {
	case METABALL_METABALL: bu_log("Sorry, strict metaballs are not yet implemented\n");
	    break;
	case METABALL_ISOPOTENTIAL:
	    for (BU_LIST_FOR(mbpt, wdb_metaballpt, &mb->metaball_ctrl_head)) {
		VSUB2(v, *p, mbpt->coord);
		a = MAGSQ(v);
		VJOIN1(*n, *n, fabs(mbpt->fldstr)*mbpt->fldstr / (SQ(a)), v);	/* f/r^4 */
	    }
	    break;
	case METABALL_BLOB:
	    for (BU_LIST_FOR(mbpt, wdb_metaballpt, &mb->metaball_ctrl_head)) {
		VSUB2(v, *p, mbpt->coord);
		a = MAGSQ(v);
		VJOIN1(*n, *n, 2.0*mbpt->sweat/SQ(mbpt->fldstr)*exp(mbpt->sweat*(1-(a/SQ(mbpt->fldstr)))) , v);
	    }
	    break;
	default: bu_log("unknown metaball method\n"); break;
    }
    VUNITIZE(*n);
}

/**
 * Makes a deep copy of a NMG model structure.
 */
struct model *
nmg_clone_model(const struct model *original)
{
    struct model           *ret;
    void *              *structArray;
    const struct nmgregion *originalRegion;
    struct bn_tol           tolerance;

    NMG_CK_MODEL(original);

    structArray = (void **)bu_calloc(original->maxindex, sizeof(void *), "nmg_clone_model() structArray");

    ret = nmg_mm();
    ret->index    = original->index;
    ret->maxindex = original->maxindex;

    structArray[ret->index] = ret;

    tolerance.magic   = BN_TOL_MAGIC;
    tolerance.dist    = 0.0005;
    tolerance.dist_sq = tolerance.dist * tolerance.dist;
    tolerance.perp    = 1e-6;
    tolerance.para    = 1 - tolerance.perp;

    for (BU_LIST_FOR(originalRegion, nmgregion, &original->r_hd)) {
	struct nmgregion *newRegion = (struct nmgregion *)structArray[originalRegion->index];

	if (newRegion == NULL) {
	    const struct shell *originalShell;

	    newRegion = nmg_construct_region(ret, originalRegion, structArray);

	    for (BU_LIST_FOR(originalShell, shell, &originalRegion->s_hd)) {
		struct shell *newShell = (struct shell *)structArray[originalShell->index];

		if (newShell == NULL)
		    newShell = nmg_construct_shell(newRegion, originalShell, structArray);

		BU_LIST_INSERT(&newRegion->s_hd, &newShell->l);
	    }

	    BU_LIST_INSERT(&ret->r_hd, &newRegion->l);
	}
    }

    bu_free(structArray, "nmg_clone_model() structArray");

    return ret;
}

struct region_pair *
add_unique_pair(struct region_pair *list, /* list to add into */
		struct region *r1,        /* first region involved */
		struct region *r2,        /* second region involved */
		double dist,              /* distance/thickness metric value */
		point_t pt)               /* location where this takes place */
{
    struct region_pair *rp, *rpair;

    /* look for it in our list */
    bu_semaphore_acquire(GED_SEM_LIST);
    for (BU_LIST_FOR (rp, region_pair, &list->l)) {

	if ((r1 == rp->r.r1 && r2 == rp->r2) || (r1 == rp->r2 && r2 == rp->r.r1)) {
	    /* we already have an entry for this region pair, we
	     * increase the counter, check the depth and update
	     * thickness maximum and entry point if need be and
	     * return.
	     */
	    rp->count++;

	    if (dist > rp->max_dist) {
		rp->max_dist = dist;
		VMOVE(rp->coord, pt);
	    }
	    rpair = rp;
	    goto found;
	}
    }
    /* didn't find it in the list.  Add it */
    BU_ALLOC(rpair, struct region_pair);
    rpair->r.r1 = r1;
    rpair->r2 = r2;
    rpair->count = 1;
    rpair->max_dist = dist;
    VMOVE(rpair->coord, pt);
    list->max_dist ++; /* really a count */

    /* insert in the list at the "nice" place */
    for (BU_LIST_FOR (rp, region_pair, &list->l)) {
	if (bu_strcmp(rp->r.r1->reg_name, r1->reg_name) <= 0)
	    break;
    }
    BU_LIST_INSERT(&rp->l, &rpair->l);
 found:
    bu_semaphore_release(GED_SEM_LIST);
    return rpair;
}

void
rt_label_vlist_verts(struct bn_vlblock *vbp, struct bu_list *src, fastf_t *mat,
        double sz, double mm2local)
{
    struct bn_vlist *vp;
    struct bu_list *vhead;
    char label[256];

    vhead = bn_vlblock_find(vbp, 255, 255, 255);	/* white */

    for (BU_LIST_FOR(vp, bn_vlist, src)) {
        register int i;
        register int nused = vp->nused;
        register int *cmd = vp->cmd;
        register point_t *pt = vp->pt;

        for (i = 0; i < nused; i++, cmd++, pt++) {
            /* default coordinates label */
            /* XXX Skip polygon markers? */
            sprintf(label, " %g, %g, %g",
                (*pt)[0]*mm2local, (*pt)[1]*mm2local, (*pt)[2]*mm2local);

            bn_vlist_3string(vhead, vbp->free_vlist_hd, label, (*pt), mat, sz);
        }
    }
}

HIDDEN struct bu_cmdhist_obj *
cho_open(ClientData UNUSED(clientData), Tcl_Interp *interp, const char *name)
{
    struct bu_cmdhist_obj *chop;

    /* check to see if command history object exists */
    for (BU_LIST_FOR(chop, bu_cmdhist_obj, &HeadCmdHistObj.l)) {
	if (BU_STR_EQUAL(name, bu_vls_addr(&chop->cho_name))) {
	    Tcl_AppendResult(interp, "ch_open: ", name,
			     " exists.\n", (char *)NULL);
	    return BU_CMDHIST_OBJ_NULL;
	}
    }

    BU_GET(chop, struct bu_cmdhist_obj);
    bu_vls_init(&chop->cho_name);
    bu_vls_strcpy(&chop->cho_name, name);
    BU_LIST_INIT(&chop->cho_head.l);
    bu_vls_init(&chop->cho_head.h_command);
    chop->cho_head.h_start.tv_sec = chop->cho_head.h_start.tv_usec =
	chop->cho_head.h_finish.tv_sec = chop->cho_head.h_finish.tv_usec = 0L;
    chop->cho_head.h_status = TCL_OK;
    chop->cho_curr = &chop->cho_head;

    BU_LIST_APPEND(&HeadCmdHistObj.l, &chop->l);
    return chop;
}

int
cho_open_tcl(ClientData clientData, Tcl_Interp *interp, int argc, const char **argv)
{
    struct bu_cmdhist_obj *chop;
    struct bu_vls vls = BU_VLS_INIT_ZERO;

    if (argc == 1) {
	/* get list of command history objects */
	for (BU_LIST_FOR(chop, bu_cmdhist_obj, &HeadCmdHistObj.l))
	    Tcl_AppendResult(interp, bu_vls_addr(&chop->cho_name), " ", (char *)NULL);

	return TCL_OK;
    }

    if (argc == 2) {
	if ((chop = cho_open(clientData, interp, argv[1])) == BU_CMDHIST_OBJ_NULL)
	    return TCL_ERROR;

	(void)Tcl_CreateCommand(interp,
				bu_vls_addr(&chop->cho_name),
				(Tcl_CmdProc *)cho_cmd,
				(ClientData)chop,
				cho_deleteProc);

	/* Return new function name as result */
	Tcl_ResetResult(interp);
	Tcl_AppendResult(interp, bu_vls_addr(&chop->cho_name), (char *)NULL);
	return TCL_OK;
    }

    bu_vls_printf(&vls, "helplib ch_open");
    Tcl_Eval(interp, bu_vls_addr(&vls));
    bu_vls_free(&vls);
    return TCL_ERROR;
}

int
nmg_eu_is_part_of_crack(const struct edgeuse *eu)
{
    struct loopuse *lu;
    struct edgeuse *eu_test;

    NMG_CK_EDGEUSE(eu);

    /* must be part of a loop to be a crack */
    if (*eu->up.magic_p != NMG_LOOPUSE_MAGIC)
	return 0;

    lu = eu->up.lu_p;
    NMG_CK_LOOPUSE(lu);

    for (BU_LIST_FOR(eu_test, edgeuse, &lu->down_hd)) {
	if (eu_test == eu)
	    continue;

	if (eu_test->vu_p->v_p == eu->eumate_p->vu_p->v_p &&
	    eu_test->eumate_p->vu_p->v_p == eu->vu_p->v_p)
	    return 1;
    }

    return 0;
}

void
rt_vlist_to_uplot(FILE *fp, const struct bu_list *vhead)
{
    register struct bn_vlist *vp;

    for (BU_LIST_FOR(vp, bn_vlist, vhead)) {
	register int i;
	register int nused = vp->nused;
	register const int *cmd = vp->cmd;
	register point_t *pt = vp->pt;

	for (i = 0; i < nused; i++, cmd++, pt++) {
	    switch (*cmd) {
		case BN_VLIST_POLY_START:
		case BN_VLIST_TRI_START:
		    break;
		case BN_VLIST_POLY_MOVE:
		case BN_VLIST_LINE_MOVE:
		case BN_VLIST_TRI_MOVE:
		    pdv_3move(fp, *pt);
		    break;
		case BN_VLIST_POLY_DRAW:
		case BN_VLIST_POLY_END:
		case BN_VLIST_LINE_DRAW:
		case BN_VLIST_TRI_DRAW:
		case BN_VLIST_TRI_END:
		    pdv_3cont(fp, *pt);
		    break;
		default:
		    bu_log("rt_vlist_to_uplot: unknown vlist cmd x%x\n",
			   *cmd);
	    }
	}
    }
}

HIDDEN int
convert_grid(int idx)
{
    struct coord_sys *cs;
    point_t tmp_pt;
    VSETALL(tmp_pt, 0.0);


    if (!g_pts[idx].cid)
	return 0;

    for (BU_LIST_FOR(cs, coord_sys, &coord_head.l)) {
	if (cs->cid != g_pts[idx].cid)
	    continue;
	break;
    }

    if (BU_LIST_IS_HEAD(&cs->l, &coord_head.l)) {
	bu_exit(1, "No coordinate system defined for grid point #%d!\n", g_pts[idx].gid);
    }

    if (convert_pt(g_pts[idx].pt, cs, tmp_pt))
	return 1;

    VMOVE(g_pts[idx].pt, tmp_pt);
    g_pts[idx].cid = 0;

    return 0;
}

void
rt_metaball_print(register const struct soltab *stp)
{
    int metaball_count = 0;
    struct rt_metaball_internal *mb;
    struct wdb_metaballpt *mbpt;

    mb = (struct rt_metaball_internal *)stp->st_specific;
    RT_METABALL_CK_MAGIC(mb);
    for (BU_LIST_FOR(mbpt, wdb_metaballpt, &mb->metaball_ctrl_head)) ++metaball_count;
    bu_log("Metaball with %d points and a threshold of %g (%s rendering)\n", metaball_count, mb->threshold, rt_metaball_lookup_type_name(mb->method));
    metaball_count = 0;
    for (BU_LIST_FOR(mbpt, wdb_metaballpt, &mb->metaball_ctrl_head))
	bu_log("\t%d: %g field strength at (%g, %g, %g) and 'goo' of %g\n", ++metaball_count, mbpt->fldstr, V3ARGS(mbpt->coord), mbpt->sweat);
    return;
}

HIDDEN void
wdb_do_inter(struct bu_list *hp)
{
    struct tokens *tok;

    for (BU_LIST_FOR(tok, tokens, hp)) {
	struct tokens *prev, *next;
	union tree *tp;

	if (tok->type != WDB_TOK_INTER)
	    continue;

	prev = BU_LIST_PREV(tokens, &tok->l);
	next = BU_LIST_NEXT(tokens, &tok->l);

	if (prev->type !=WDB_TOK_TREE || next->type != WDB_TOK_TREE)
	    continue;

	/* this is an eligible intersection operation */
	BU_ALLOC(tp, union tree);
	RT_TREE_INIT(tp);
	tp->tr_b.tb_op = OP_INTERSECT;
	tp->tr_b.tb_regionp = (struct region *)NULL;
	tp->tr_b.tb_left = prev->tp;
	tp->tr_b.tb_right = next->tp;
	BU_LIST_DEQUEUE(&tok->l);
	bu_free((char *)tok, "tok");
	BU_LIST_DEQUEUE(&prev->l);
	bu_free((char *)prev, "prev");
	next->tp = tp;
	tok = next;
    }
}

HIDDEN int
wdb_do_paren(struct bu_list *hp)
{
    struct tokens *tok;

    for (BU_LIST_FOR(tok, tokens, hp)) {
	struct tokens *prev, *next;

	if (tok->type != WDB_TOK_TREE)
	    continue;

	prev = BU_LIST_PREV(tokens, &tok->l);
	next = BU_LIST_NEXT(tokens, &tok->l);

	if (prev->type !=WDB_TOK_LPAREN || next->type != WDB_TOK_RPAREN)
	    continue;

	/* this is an eligible operand surrounded by parens */
	BU_LIST_DEQUEUE(&next->l);
	bu_free((char *)next, "next");
	BU_LIST_DEQUEUE(&prev->l);
	bu_free((char *)prev, "prev");
    }

    if (hp->forw == hp->back && hp->forw != hp)
	return 1;	/* done */
    else if (BU_LIST_IS_EMPTY(hp))
	return -1;	/* empty tree!!!! */
    else
	return 0;	/* more to do */

}

static void
Write_euclid_face(const struct loopuse *lu, const int facet_type, const int regionid, const int face_number)
{
    struct faceuse *fu;
    struct edgeuse *eu;
    plane_t plane;
    int vertex_count = 0;

    NMG_CK_LOOPUSE(lu);

    if (verbose)
	bu_log("Write_euclid_face: lu=%p, facet_type=%d, regionid=%d, face_number=%d\n",
	       (void *)lu, facet_type, regionid, face_number);

    if (BU_LIST_FIRST_MAGIC(&lu->down_hd) != NMG_EDGEUSE_MAGIC)
	return;

    if (*lu->up.magic_p != NMG_FACEUSE_MAGIC)
	return;

    for (BU_LIST_FOR(eu, edgeuse, &lu->down_hd))
	vertex_count++;

    fprintf(fp_out, "%10d%3d     0.    1%5d", regionid, facet_type, vertex_count);

    vertex_count = 0;
    for (BU_LIST_FOR(eu, edgeuse, &lu->down_hd))
    {
	struct vertex *v;
	int i;

	NMG_CK_EDGEUSE(eu);
	v = eu->vu_p->v_p;
	NMG_CK_VERTEX(v);
/*		fprintf(fp_out, "%10d%8f%8f%8f", ++vertex_count, V3ARGS(v->vg_p->coord)); */
	vertex_count++;
	fprintf(fp_out, "%10d", vertex_count);

	for (i=X; i<=Z; i++)
	    fastf_print(fp_out, 8, v->vg_p->coord[i]);
    }

    fu = lu->up.fu_p;
    NMG_CK_FACEUSE(fu);
    NMG_GET_FU_PLANE(plane, fu);
    fprintf(fp_out, "%10d%15.5f%15.5f%15.5f%15.5f", face_number, V4ARGS(plane));
}

/**
 * Intersect all the edges in fu1 that don't lie on any of the faces
 * of shell s2 with s2, i.e. "interior" edges, where the endpoints lie
 * on s2, but the edge is not shared with a face of s2.  Such edges
 * wouldn't have been processed by the NEWLINE version of
 * nmg_isect_two_generic_faces(), so intersections need to be looked
 * for here.  Fortunately, it's easy to reject everything except edges
 * that need processing using only the topology structures.
 *
 * The "_int" at the end of the name is to signify that this routine
 * does only "interior" edges, and is not a general face/shell
 * intersector.
 */
void
nmg_isect_face3p_shell_int(struct nmg_inter_struct *is, struct faceuse *fu1, struct shell *s2)
{
    struct shell *s1;
    struct loopuse *lu1;
    struct edgeuse *eu1;

    NMG_CK_INTER_STRUCT(is);
    NMG_CK_FACEUSE(fu1);
    NMG_CK_SHELL(s2);
    s1 = fu1->s_p;
    NMG_CK_SHELL(s1);

    if (RTG.NMG_debug & DEBUG_POLYSECT)
	bu_log("nmg_isect_face3p_shell_int(, fu1=x%x, s2=x%x) START\n", fu1, s2);

    for (BU_LIST_FOR (lu1, loopuse, &fu1->lu_hd)) {
	NMG_CK_LOOPUSE(lu1);
	if (BU_LIST_FIRST_MAGIC(&lu1->down_hd) == NMG_VERTEXUSE_MAGIC)
	    continue;
	for (BU_LIST_FOR (eu1, edgeuse, &lu1->down_hd)) {
	    struct edgeuse *eu2;

	    eu2 = nmg_find_matching_eu_in_s(eu1, s2);
	    if (eu2) {
		bu_log("nmg_isect_face3p_shell_int() eu1=x%x, e1=x%x, eu2=x%x, e2=x%x (nothing to do)\n", eu1, eu1->e_p, eu2, eu2->e_p);
		/* Whether the edgeuse is in a face, or a wire
		 * edgeuse, the other guys will isect it.
		 */
		continue;
	    }
	    /* vu2a and vu2b are in shell s2, but there is no edge
	     * running between them in shell s2.  Create a line of
	     * intersection, and go to it!.
	     */
	    bu_log("nmg_isect_face3p_shell_int(, s2=x%x) eu1=x%x, no eu2\n", s2, eu1);
/* XXX eso no existe todavia */
/* 	    nmg_isect_edge3p_shell(is, eu1, s2); */
	}
    }

    if (RTG.NMG_debug & DEBUG_POLYSECT)
	bu_log("nmg_isect_face3p_shell_int(, fu1=x%x, s2=x%x) END\n", fu1, s2);
}

/**
 * Mesh every edge in shell 1 with every edge in shell 2.
 * The return is the number of edges meshed.
 *
 * Does not use nmg_mesh_face_shell() to keep face/self meshing
 * to the absolute minimum necessary.
 */
int
nmg_mesh_shell_shell(struct shell *s1, struct shell *s2, const struct bn_tol *tol)
{
    struct faceuse *fu1;
    struct faceuse *fu2;
    int count = 0;

    NMG_CK_SHELL(s1);
    NMG_CK_SHELL(s2);
    BN_CK_TOL(tol);

    nmg_region_v_unique(s1->r_p, tol);
    nmg_region_v_unique(s2->r_p, tol);

    /* First, mesh all faces of shell 2 with themselves */
    for (BU_LIST_FOR(fu2, faceuse, &s2->fu_hd)) {
	NMG_CK_FACEUSE(fu2);
	count += nmg_mesh_two_faces(fu2, fu2, tol);
    }

    /* Visit every face in shell 1 */
    for (BU_LIST_FOR(fu1, faceuse, &s1->fu_hd)) {
	NMG_CK_FACEUSE(fu1);

	/* First, mesh each face in shell 1 with itself */
	count += nmg_mesh_two_faces(fu1, fu1, tol);

	/* Visit every face in shell 2 */
	for (BU_LIST_FOR(fu2, faceuse, &s2->fu_hd)) {
	    NMG_CK_FACEUSE(fu2);
	    count += nmg_mesh_two_faces(fu1, fu2, tol);
	}
    }

    /* XXX What about wire edges in the shell? */

    /* Visit every wire loop in shell 1 */

    /* Visit every wire edge in shell 1 */

    return count;
}

HIDDEN void
get_cbar(void)
{
    int eid, pid;
    int g1, g2;
    point_t pt1, pt2;
    fastf_t radius;
    vect_t height;
    struct pbar *pb;
    char cbar_name[NAMESIZE+1];

    eid = atoi( curr_rec[1] );

    pid = atoi( curr_rec[2] );
    if ( !pid )
    {
	if ( bar_def_pid )
	    pid = bar_def_pid;
	else
	    pid = eid;
    }

    g1 = atoi( curr_rec[3] );

    g2 = atoi( curr_rec[4] );

    get_grid( g1, pt1 );
    get_grid( g2, pt2 );

    for ( BU_LIST_FOR( pb, pbar, &pbar_head.l ) )
    {
	if ( pb->pid == pid )
	    break;
    }

    if ( BU_LIST_IS_HEAD( &pb->l, &pbar_head.l ) )
    {
	log_line( "Non-existent PID referenced in CBAR" );
	return;
    }

    VSCALE( pt1, pt1, conv[units] );
    VSCALE( pt2, pt2, conv[units] );

    radius = sqrt( pb->area/bn_pi );
    radius = radius * conv[units];

    VSUB2( height, pt2, pt1 );

    sprintf( cbar_name, "cbar.%d", eid );
    mk_rcc( fpout, cbar_name, pt1, height, radius );

    mk_addmember( cbar_name, &pb->head.l, NULL, WMOP_UNION );
}

void show_sites (struct bu_list *sl)
{
    struct site *sp;

    BU_CK_LIST_HEAD(sl);

    for (BU_LIST_FOR(sp, site, sl)) {
	bu_log("I got a site (%g, %g, %g)\n",
	       sp->s_x, sp->s_y, sp->s_z);
    }
}