src/hd/rhoptimize.c

/* Copyright (c) 1998 Silicon Graphics, Inc. */

#ifndef lint
static char SCCSid[] = "$SunId$ SGI";
#endif

/*
 * Optimize holodeck for quick access.
 *
 *      11/4/98         Greg Ward Larson
 */

#include "holo.h"

#ifndef BKBSIZE
#define BKBSIZE         256             /* beam clump size (kilobytes) */
#endif

#define flgop(p,i,op)           ((p)[(i)>>5] op (1L<<((i)&0x1f)))
#define isset(p,i)              flgop(p,i,&)
#define setfl(p,i)              flgop(p,i,|=)
#define clrfl(p,i)              flgop(p,i,&=~)

char    *progname;

extern long     rhinitcopy();


main(argc, argv)
int     argc;
char    *argv[];
{
        char    nambuf[128];
        char    *inpname, *outname;
        int     hdfd[2];
        long    nextipos, lastopos, thisopos;

        progname = argv[0];
        if (argc < 2 | argc > 3) {
                fprintf(stderr, "Usage: %s input.hdk [output.hdk]\n", progname);
                exit(1);
        }
        inpname = argv[1];
        if (argc == 3)                  /* use given output file */
                outname = argv[2];
        else {                          /* else create temporary file */
                strcpy(nambuf, inpname);
                if ((outname = strrchr(nambuf, '/')) != NULL)
                        outname++;
                else
                        outname = nambuf;
                sprintf(outname, "rho%d.hdk", getpid());
                outname = nambuf;
        }
                                        /* copy holodeck file header */
        nextipos = rhinitcopy(hdfd, inpname, outname);
        lastopos = 0L;                  /* copy sections one by one */
        while (nextipos != 0L) {
                                        /* set input position; get next */
                lseek(hdfd[0], nextipos, 0);
                read(hdfd[0], (char *)&nextipos, sizeof(nextipos));
                                        /* get output position; set last */
                thisopos = lseek(hdfd[1], 0L, 2);
                if (lastopos > 0L) {
                        lseek(hdfd[1], lastopos, 0);
                        write(hdfd[1], (char *)&thisopos, sizeof(thisopos));
                        lseek(hdfd[1], 0L, 2);
                }
                lastopos = thisopos;
                thisopos = 0L;          /* write place holder */
                write(hdfd[1], (char *)&thisopos, sizeof(thisopos));
                                        /* copy holodeck section */
                copysect(hdfd[0], hdfd[1]);
        }
                                        /* clean up */
        close(hdfd[0]);
        close(hdfd[1]);
        if (argc == 2 && rename(outname, inpname) < 0) {
                sprintf(errmsg, "cannot rename \"%s\" to \"%s\"",
                                outname, inpname);
                error(SYSTEM, errmsg);
        }
        exit(0);
}


long
rhinitcopy(hfd, infn, outfn)    /* open files and copy header */
int     hfd[2];                 /* returned file descriptors */
char    *infn, *outfn;
{
        FILE    *infp, *outfp;
        long    ifpos;
                                        /* open files for i/o */
        if ((infp = fopen(infn, "r")) == NULL) {
                sprintf(errmsg, "cannot open \"%s\" for reading", infn);
                error(SYSTEM, errmsg);
        }
        if ((outfp = fopen(outfn, "w+")) == NULL) {
                sprintf(errmsg, "cannot open \"%s\" for writing", outfn);
                error(SYSTEM, errmsg);
        }
                                        /* copy and verify header */
        if (checkheader(infp, HOLOFMT, outfp) < 0 ||
                        getw(infp) != HOLOMAGIC)
                error(USER, "input not in holodeck format");
        fputformat(HOLOFMT, outfp);
        fputc('\n', outfp);
        putw(HOLOMAGIC, outfp);
                                        /* get descriptors and free stdio */
        if ((hfd[0] = dup(fileno(infp))) < 0 ||
                        (hfd[1] = dup(fileno(outfp))) < 0)
                error(SYSTEM, "dup call failed in rhinitcopy");
        ifpos = ftell(infp);
        fclose(infp);
        if (fclose(outfp) == EOF)
                error(SYSTEM, "file flushing error in rhinitcopy");
                                        /* we flush everything manually */
        hdcachesize = 0;
                                        /* return input position */
        return(ifpos);
}


gcshifti(gc, ia, di, hp)        /* shift cell row or column */
register GCOORD *gc;
int     ia, di;
register HOLO   *hp;
{
        int     nw;

        if (di > 0) {
                if (++gc->i[ia] >= hp->grid[((gc->w>>1)+1+ia)%3]) {
                        nw = ((gc->w&~1) + (ia<<1) + 3) % 6;
                        gc->i[ia] = gc->i[1-ia];
                        gc->i[1-ia] = gc->w&1 ? hp->grid[((nw>>1)+2-ia)%3]-1 : 0;
                        gc->w = nw;
                }
        } else if (di < 0) {
                if (--gc->i[ia] < 0) {
                        nw = ((gc->w&~1) + (ia<<1) + 2) % 6;
                        gc->i[ia] = gc->i[1-ia];
                        gc->i[1-ia] = gc->w&1 ? hp->grid[((nw>>1)+2-ia)%3]-1 : 0;
                        gc->w = nw;
                }
        }
}


mkneighgrid(ng, hp, gc)         /* compute neighborhood for grid cell */
GCOORD  ng[3*3];
HOLO    *hp;
GCOORD  *gc;
{
        GCOORD  gci0;
        register int    i, j;

        for (i = 3; i--; ) {
                copystruct(&gci0, gc);
                gcshifti(&gci0, 0, i-1, hp);
                for (j = 3; j--; ) {
                        copystruct(ng+(3*i+j), &gci0);
                        gcshifti(ng+(3*i+j), gci0.w==gc->w, j-1, hp);
                }
        }
}


int     bneighlist[9*9-1];
int     bneighrem;

#define nextneigh()     (bneighrem<=0 ? 0 : bneighlist[--bneighrem])

int
firstneigh(hp, b)               /* initialize neighbor list and return first */
HOLO    *hp;
int     b;
{
        GCOORD  wg0[9], wg1[9], bgc[2];
        int     i, j;

        hdbcoord(bgc, hp, b);
        mkneighgrid(wg0, hp, bgc);
        mkneighgrid(wg1, hp, bgc+1);
        bneighrem = 0;
        for (i = 9; i--; )
                for (j = 9; j--; ) {
                        if (i == 4 & j == 4)
                                continue;
                        if (wg0[i].w == wg1[j].w)
                                continue;
                        copystruct(bgc, wg0+i);
                        copystruct(bgc+1, wg1+j);
                        bneighlist[bneighrem++] = hdbindex(hp, bgc);
#ifdef DEBUG
                        if (bneighlist[bneighrem-1] <= 0)
                                error(CONSISTENCY, "bad beam in firstneigh");
#endif
                }
        return(nextneigh());
}


BEAMI   *beamdir;

int
bpcmp(b1p, b2p)                 /* compare beam positions on disk */
int     *b1p, *b2p;
{
        register long   pdif = beamdir[*b1p].fo - beamdir[*b2p].fo;

        if (pdif > 0) return(1);
        if (pdif < 0) return(-1);
        return(0);
}


copysect(ifd, ofd)              /* copy holodeck section from ifd to ofd */
int     ifd, ofd;
{
        static short    primes[] = {9431,6803,4177,2659,1609,887,587,251,47,1};
        register HOLO   *hinp;
        HOLO    *hout;
        register BEAM   *bp;
        unsigned int4   *bflags;
        int     *bqueue;
        int     bqlen;
        int4    bqtotal;
        int     bc, bci, bqc, myprime;
        register int    i;
                                        /* load input section directory */
        hinp = hdinit(ifd, NULL);
                                        /* create output section directory */
        hout = hdinit(ofd, (HDGRID *)hinp);
                                        /* allocate beam queue */
        bqueue = (int *)malloc(nbeams(hinp)*sizeof(int));
        bflags = (int4 *)calloc((nbeams(hinp)>>3)+1, sizeof(int4));
        if (bqueue == NULL | bflags == NULL)
                error(SYSTEM, "out of memory in copysect");
                                        /* mark empty beams as done */
        for (i = nbeams(hinp); i-- > 0; )
                if (!hinp->bi[i].nrd)
                        setfl(bflags, i);
                                        /* pick a good prime step size */
        for (i = 0; primes[i]<<5 >= nbeams(hinp); i++)
                ;
        while ((myprime = primes[i++]) > 1)
                if (nbeams(hinp) % myprime)
                        break;
                                        /* add each input beam and neighbors */
        for (bc = bci = nbeams(hinp); bc > 0; bc--,
                        bci += bci>myprime ? -myprime : nbeams(hinp)-myprime) {
                if (isset(bflags, bci))
                        continue;
                bqueue[0] = bci;                /* initialize queue */
                bqlen = 1;
                bqtotal = bnrays(hinp, bci);
                setfl(bflags, bci);
                                                /* run through growing queue */
                for (bqc = 0; bqc < bqlen; bqc++) {
                                                /* add neighbors until full */
                        for (i = firstneigh(hinp,bqueue[bqc]); i > 0;
                                        i = nextneigh()) {
                                if (isset(bflags, i))   /* done already? */
                                        continue;
                                bqueue[bqlen++] = i;    /* add it */
                                bqtotal += bnrays(hinp, i);
                                setfl(bflags, i);
                                if (bqtotal >= BKBSIZE*1024/sizeof(RAYVAL))
                                        break;          /* queue full */
                        }
                        if (i > 0)
                                break;
                }
                beamdir = hinp->bi;             /* sort queue */
                qsort((char *)bqueue, bqlen, sizeof(*bqueue), bpcmp);
                                                /* transfer each beam */
                for (i = 0; i < bqlen; i++) {
                        bp = hdgetbeam(hinp, bqueue[i]);
                        bcopy((char *)hdbray(bp),
                                (char *)hdnewrays(hout,bqueue[i],bp->nrm),
                                        bp->nrm*sizeof(RAYVAL));
                        hdfreebeam(hinp, bqueue[i]);
                }
                hdfreebeam(hout, 0);            /* flush output block */
#ifdef DEBUG
                hdsync(hout, 0);
#endif
        }
                                        /* we're done -- clean up */
        free((char *)bqueue);
        free((char *)bflags);
        hddone(hinp);
        hddone(hout);
}


eputs(s)                        /* put error message to stderr */
register char  *s;
{
        static int  midline = 0;

        if (!*s)
                return;
        if (!midline++) {       /* prepend line with program name */
                fputs(progname, stderr);
                fputs(": ", stderr);
        }
        fputs(s, stderr);
        if (s[strlen(s)-1] == '\n') {
                fflush(stderr);
                midline = 0;
        }
}


quit(code)                      /* exit the program gracefully */
int     code;
{
        hdsync(NULL, 1);        /* write out any buffered data */
        exit(code);
}
Revision:	3.2
Committed:	Thu Nov 5 14:43:14 1998 UTC (26 years, 11 months ago) by gwlarson
Content type:	text/plain
Branch:	MAIN
Changes since 3.1:	+63 -28 lines
Log Message:	build queue then read bug fixes in neighbor algorithm
#	User	Rev	Content
1	gwlarson	3.1	/* Copyright (c) 1998 Silicon Graphics, Inc. */
2
3			#ifndef lint
4			static char SCCSid[] = "$SunId$ SGI";
5			#endif
6
7			/*
8			* Optimize holodeck for quick access.
9			*
10			* 11/4/98 Greg Ward Larson
11			*/
12
13			#include "holo.h"
14
15			#ifndef BKBSIZE
16	gwlarson	3.2	#define BKBSIZE 256 /* beam clump size (kilobytes) */
17	gwlarson	3.1	#endif
18
19	gwlarson	3.2	#define flgop(p,i,op) ((p)[(i)>>5] op (1L<<((i)&0x1f)))
20			#define isset(p,i) flgop(p,i,&)
21			#define setfl(p,i) flgop(p,i,\|=)
22			#define clrfl(p,i) flgop(p,i,&=~)
23	gwlarson	3.1
24			char *progname;
25
26			extern long rhinitcopy();
27
28
29			main(argc, argv)
30			int argc;
31			char *argv[];
32			{
33			char nambuf[128];
34			char inpname, outname;
35			int hdfd[2];
36			long nextipos, lastopos, thisopos;
37
38			progname = argv[0];
39			if (argc < 2 \| argc > 3) {
40			fprintf(stderr, "Usage: %s input.hdk [output.hdk]\n", progname);
41			exit(1);
42			}
43			inpname = argv[1];
44			if (argc == 3) /* use given output file */
45			outname = argv[2];
46			else { /* else create temporary file */
47			strcpy(nambuf, inpname);
48			if ((outname = strrchr(nambuf, '/')) != NULL)
49			outname++;
50			else
51			outname = nambuf;
52			sprintf(outname, "rho%d.hdk", getpid());
53			outname = nambuf;
54			}
55			/* copy holodeck file header */
56			nextipos = rhinitcopy(hdfd, inpname, outname);
57			lastopos = 0L; /* copy sections one by one */
58			while (nextipos != 0L) {
59			/* set input position; get next */
60			lseek(hdfd[0], nextipos, 0);
61			read(hdfd[0], (char *)&nextipos, sizeof(nextipos));
62			/* get output position; set last */
63			thisopos = lseek(hdfd[1], 0L, 2);
64			if (lastopos > 0L) {
65			lseek(hdfd[1], lastopos, 0);
66			write(hdfd[1], (char *)&thisopos, sizeof(thisopos));
67			lseek(hdfd[1], 0L, 2);
68			}
69			lastopos = thisopos;
70			thisopos = 0L; /* write place holder */
71			write(hdfd[1], (char *)&thisopos, sizeof(thisopos));
72			/* copy holodeck section */
73			copysect(hdfd[0], hdfd[1]);
74			}
75			/* clean up */
76			close(hdfd[0]);
77			close(hdfd[1]);
78			if (argc == 2 && rename(outname, inpname) < 0) {
79			sprintf(errmsg, "cannot rename \"%s\" to \"%s\"",
80			outname, inpname);
81			error(SYSTEM, errmsg);
82			}
83			exit(0);
84			}
85
86
87			long
88			rhinitcopy(hfd, infn, outfn) /* open files and copy header */
89			int hfd[2]; /* returned file descriptors */
90			char infn, outfn;
91			{
92			FILE infp, outfp;
93			long ifpos;
94			/* open files for i/o */
95			if ((infp = fopen(infn, "r")) == NULL) {
96			sprintf(errmsg, "cannot open \"%s\" for reading", infn);
97			error(SYSTEM, errmsg);
98			}
99			if ((outfp = fopen(outfn, "w+")) == NULL) {
100			sprintf(errmsg, "cannot open \"%s\" for writing", outfn);
101			error(SYSTEM, errmsg);
102			}
103			/* copy and verify header */
104			if (checkheader(infp, HOLOFMT, outfp) < 0 \|\|
105			getw(infp) != HOLOMAGIC)
106			error(USER, "input not in holodeck format");
107			fputformat(HOLOFMT, outfp);
108			fputc('\n', outfp);
109			putw(HOLOMAGIC, outfp);
110			/* get descriptors and free stdio */
111			if ((hfd[0] = dup(fileno(infp))) < 0 \|\|
112			(hfd[1] = dup(fileno(outfp))) < 0)
113			error(SYSTEM, "dup call failed in rhinitcopy");
114			ifpos = ftell(infp);
115			fclose(infp);
116			if (fclose(outfp) == EOF)
117			error(SYSTEM, "file flushing error in rhinitcopy");
118			/* we flush everything manually */
119			hdcachesize = 0;
120			/* return input position */
121			return(ifpos);
122			}
123
124
125			gcshifti(gc, ia, di, hp) /* shift cell row or column */
126			register GCOORD *gc;
127			int ia, di;
128			register HOLO *hp;
129			{
130			int nw;
131
132			if (di > 0) {
133			if (++gc->i[ia] >= hp->grid[((gc->w>>1)+1+ia)%3]) {
134			nw = ((gc->w&~1) + (ia<<1) + 3) % 6;
135			gc->i[ia] = gc->i[1-ia];
136			gc->i[1-ia] = gc->w&1 ? hp->grid[((nw>>1)+2-ia)%3]-1 : 0;
137			gc->w = nw;
138			}
139			} else if (di < 0) {
140			if (--gc->i[ia] < 0) {
141			nw = ((gc->w&~1) + (ia<<1) + 2) % 6;
142			gc->i[ia] = gc->i[1-ia];
143			gc->i[1-ia] = gc->w&1 ? hp->grid[((nw>>1)+2-ia)%3]-1 : 0;
144			gc->w = nw;
145			}
146			}
147			}
148
149
150			mkneighgrid(ng, hp, gc) /* compute neighborhood for grid cell */
151			GCOORD ng[3*3];
152			HOLO *hp;
153			GCOORD *gc;
154			{
155			GCOORD gci0;
156	gwlarson	3.2	register int i, j;
157	gwlarson	3.1
158			for (i = 3; i--; ) {
159			copystruct(&gci0, gc);
160			gcshifti(&gci0, 0, i-1, hp);
161			for (j = 3; j--; ) {
162			copystruct(ng+(3*i+j), &gci0);
163	gwlarson	3.2	gcshifti(ng+(3*i+j), gci0.w==gc->w, j-1, hp);
164	gwlarson	3.1	}
165			}
166			}
167
168
169			int bneighlist[9*9-1];
170			int bneighrem;
171
172			#define nextneigh() (bneighrem<=0 ? 0 : bneighlist[--bneighrem])
173
174			int
175			firstneigh(hp, b) /* initialize neighbor list and return first */
176			HOLO *hp;
177			int b;
178			{
179			GCOORD wg0[9], wg1[9], bgc[2];
180			int i, j;
181
182			hdbcoord(bgc, hp, b);
183			mkneighgrid(wg0, hp, bgc);
184			mkneighgrid(wg1, hp, bgc+1);
185			bneighrem = 0;
186			for (i = 9; i--; )
187			for (j = 9; j--; ) {
188			if (i == 4 & j == 4)
189			continue;
190	gwlarson	3.2	if (wg0[i].w == wg1[j].w)
191			continue;
192	gwlarson	3.1	copystruct(bgc, wg0+i);
193			copystruct(bgc+1, wg1+j);
194			bneighlist[bneighrem++] = hdbindex(hp, bgc);
195	gwlarson	3.2	#ifdef DEBUG
196			if (bneighlist[bneighrem-1] <= 0)
197			error(CONSISTENCY, "bad beam in firstneigh");
198			#endif
199	gwlarson	3.1	}
200			return(nextneigh());
201			}
202
203
204	gwlarson	3.2	BEAMI *beamdir;
205
206			int
207			bpcmp(b1p, b2p) /* compare beam positions on disk */
208			int b1p, b2p;
209			{
210			register long pdif = beamdir[b1p].fo - beamdir[b2p].fo;
211
212			if (pdif > 0) return(1);
213			if (pdif < 0) return(-1);
214			return(0);
215			}
216
217
218	gwlarson	3.1	copysect(ifd, ofd) /* copy holodeck section from ifd to ofd */
219			int ifd, ofd;
220			{
221			static short primes[] = {9431,6803,4177,2659,1609,887,587,251,47,1};
222	gwlarson	3.2	register HOLO *hinp;
223			HOLO *hout;
224	gwlarson	3.1	register BEAM *bp;
225	gwlarson	3.2	unsigned int4 *bflags;
226	gwlarson	3.1	int *bqueue;
227			int bqlen;
228			int4 bqtotal;
229	gwlarson	3.2	int bc, bci, bqc, myprime;
230	gwlarson	3.1	register int i;
231			/* load input section directory */
232			hinp = hdinit(ifd, NULL);
233			/* create output section directory */
234			hout = hdinit(ofd, (HDGRID *)hinp);
235			/* allocate beam queue */
236	gwlarson	3.2	bqueue = (int )malloc(nbeams(hinp)sizeof(int));
237			bflags = (int4 *)calloc((nbeams(hinp)>>3)+1, sizeof(int4));
238			if (bqueue == NULL \| bflags == NULL)
239	gwlarson	3.1	error(SYSTEM, "out of memory in copysect");
240	gwlarson	3.2	/* mark empty beams as done */
241			for (i = nbeams(hinp); i-- > 0; )
242			if (!hinp->bi[i].nrd)
243			setfl(bflags, i);
244	gwlarson	3.1	/* pick a good prime step size */
245			for (i = 0; primes[i]<<5 >= nbeams(hinp); i++)
246			;
247			while ((myprime = primes[i++]) > 1)
248			if (nbeams(hinp) % myprime)
249			break;
250			/* add each input beam and neighbors */
251			for (bc = bci = nbeams(hinp); bc > 0; bc--,
252			bci += bci>myprime ? -myprime : nbeams(hinp)-myprime) {
253	gwlarson	3.2	if (isset(bflags, bci))
254	gwlarson	3.1	continue;
255			bqueue[0] = bci; /* initialize queue */
256			bqlen = 1;
257			bqtotal = bnrays(hinp, bci);
258	gwlarson	3.2	setfl(bflags, bci);
259	gwlarson	3.1	/* run through growing queue */
260			for (bqc = 0; bqc < bqlen; bqc++) {
261	gwlarson	3.2	/* add neighbors until full */
262			for (i = firstneigh(hinp,bqueue[bqc]); i > 0;
263			i = nextneigh()) {
264			if (isset(bflags, i)) /* done already? */
265	gwlarson	3.1	continue;
266	gwlarson	3.2	bqueue[bqlen++] = i; /* add it */
267			bqtotal += bnrays(hinp, i);
268			setfl(bflags, i);
269	gwlarson	3.1	if (bqtotal >= BKBSIZE*1024/sizeof(RAYVAL))
270			break; /* queue full */
271			}
272	gwlarson	3.2	if (i > 0)
273			break;
274	gwlarson	3.1	}
275	gwlarson	3.2	beamdir = hinp->bi; /* sort queue */
276			qsort((char )bqueue, bqlen, sizeof(bqueue), bpcmp);
277			/* transfer each beam */
278			for (i = 0; i < bqlen; i++) {
279			bp = hdgetbeam(hinp, bqueue[i]);
280			bcopy((char *)hdbray(bp),
281			(char *)hdnewrays(hout,bqueue[i],bp->nrm),
282			bp->nrm*sizeof(RAYVAL));
283			hdfreebeam(hinp, bqueue[i]);
284			}
285	gwlarson	3.1	hdfreebeam(hout, 0); /* flush output block */
286	gwlarson	3.2	#ifdef DEBUG
287			hdsync(hout, 0);
288			#endif
289	gwlarson	3.1	}
290			/* we're done -- clean up */
291			free((char *)bqueue);
292	gwlarson	3.2	free((char *)bflags);
293	gwlarson	3.1	hddone(hinp);
294			hddone(hout);
295			}
296
297
298			eputs(s) /* put error message to stderr */
299			register char *s;
300			{
301			static int midline = 0;
302
303			if (!*s)
304			return;
305			if (!midline++) { /* prepend line with program name */
306			fputs(progname, stderr);
307			fputs(": ", stderr);
308			}
309			fputs(s, stderr);
310			if (s[strlen(s)-1] == '\n') {
311			fflush(stderr);
312			midline = 0;
313			}
314			}
315
316
317			quit(code) /* exit the program gracefully */
318			int code;
319			{
320			hdsync(NULL, 1); /* write out any buffered data */
321			exit(code);
322			}