ximage.c

地震波正演和显示模块

#include "../include/par.h"
#include "../include/xplot.h"
#include "../lib/getpars.c"
#include "../lib/atopkge.c"
#include "../lib/sort.c"
#include "../lib/colormaps.c"
#include "../lib/rubberbox.c"
#include "../lib/image.c"
#include "../lib/legendbox.c"
#include "../lib/axesbox.c"
#include "../lib/drawcurve.c"


/*********************** self documentation **********************/
char *sdoc[] = {
"									",
" XIMAGE - X IMAGE plot of a uniformly-sampled function f(x1,x2)       	",
"									",
" ximage n1= [optional parameters] <binaryfile				",
"									",
" X Functionality:							",
" Button 1	Zoom with rubberband box				",
" Button 2	Show mouse (x1,x2) coordinates while pressed		",
" q or Q key	Quit							",
" s key		Save current mouse (x1,x2) location to file		",
" p or P key	Plot current window with pswigb (only from disk files)	",
"									",
" ... change colormap interactively					",
" r	     install next RGB - colormap				",
" R	     install previous RGB - colormap				",
" h	     install next HSV - colormap				",
" H	     install previous HSV - colormap				",
" H	     install previous HSV - colormap				",
" (Move mouse cursor out and back into window for r,R,h,H to take effect)",
" 									",
" Required Parameters:							",
" n1			 number of samples in 1st (fast) dimension	",
"									",
" Optional Parameters:							",
" d1=1.0		 sampling interval in 1st dimension		",
" f1=0.0		 first sample in 1st dimension			",
" n2=all		 number of samples in 2nd (slow) dimension	",
" d2=1.0		 sampling interval in 2nd dimension		",
" f2=0.0		 first sample in 2nd dimension			",
" mpicks=/dev/tty	 file to save mouse picks in			",
" perc=100.0		 percentile used to determine clip		",
" clip=(perc percentile) clip used to determine bclip and wclip		",
" bperc=perc		 percentile for determining black clip value	",
" wperc=100.0-perc	 percentile for determining white clip value	",
" bclip=clip		 data values outside of [bclip,wclip] are clipped",
" wclip=-clip		 data values outside of [bclip,wclip] are clipped",
" cmap=hsv\'n\' or rgb\'m\'	\'n\' is a number from 0 to 13		",
"				\'m\' is a number from 0 to 11		",
"				cmap=rgb0 is equal to cmap=gray		",
"				cmap=hsv1 is equal to cmap=hue		",
"				(compatibility to older versions)	",
" legend=0	      =1 display the color scale			",
" units=		unit label for legend				",
" legendfont=times_roman10    font name for title			",
" verbose=1		 =1 for info printed on stderr (0 for no info)	",
" xbox=50		 x in pixels of upper left corner of window	",
" ybox=50		 y in pixels of upper left corner of window	",
" wbox=550		 width in pixels of window			",
" hbox=700		 height in pixels of window			",
" lwidth=16		 colorscale (legend) width in pixels		",
" lheight=hbox/3	 colorscale (legend) height in pixels		",
" lx=3			 colorscale (legend) x-position in pixels	",
" ly=(hbox-lheight)/3    colorscale (legend) y-position in pixels	",
" x1beg=x1min		 value at which axis 1 begins			",
" x1end=x1max		 value at which axis 1 ends			",
" d1num=0.0		 numbered tic interval on axis 1 (0.0 for automatic)",
" f1num=x1min		 first numbered tic on axis 1 (used if d1num not 0.0)",
" n1tic=1		 number of tics per numbered tic on axis 1	",
" grid1=none		 grid lines on axis 1 - none, dot, dash, or solid",
" label1=		 label on axis 1				",
" x2beg=x2min		 value at which axis 2 begins			",
" x2end=x2max		 value at which axis 2 ends			",
" d2num=0.0		 numbered tic interval on axis 2 (0.0 for automatic)",
" f2num=x2min		 first numbered tic on axis 2 (used if d2num not 0.0)",
" n2tic=1		 number of tics per numbered tic on axis 2	",
" grid2=none		 grid lines on axis 2 - none, dot, dash, or solid",
" label2=		 label on axis 2				",
" labelfont=Erg14	 font name for axes labels			",
" title=		 title of plot					",
" titlefont=Rom22	 font name for title				",
" windowtitle=ximage	 title on window				",
" labelcolor=blue	 color for axes labels				",
" titlecolor=red	 color for title				",
" gridcolor=blue	 color for grid lines				",
" style=seismic		 normal (axis 1 horizontal, axis 2 vertical) or ",
"			 seismic (axis 1 vertical, axis 2 horizontal)	",
" blank=0		 This indicates what portion of the lower range ",
"			 to blank out (make the background color).  The ",
"			 value should range from 0 to 1.		",
" curve=curve1,curve2,...  file(s) containing points to draw curve(s)   ",
" npair=n1,n2,n2,...            number(s) of pairs in each file         ",
" curvecolor=color1,color2,...  color(s) for curve(s)                   ",
" blockinterp=0       whether to use block interpolation (0=no, 1=yes)  ",
"									",
"									",
" NOTES:								",
" Currently, the curve file must have the number of points specified    ",
" at the top of the file.  This file is an ascii file with the points   ",
" specified as (x1,x2) pairs, one pair to a line.  A \"vector\" of curve",
" files and curve colors may be specified as curvefile=file1,file2,etc. ",
" and curvecolor=color1,color2,etc, and the number of pairs of values   ",
" in each file as npair=npair1,npair2,... .                             ",
"									",
NULL};

/**************** end self doc ********************************/

/* functions defined and used internally */
static void zoomBox (int x, int y, int w, int h, 
	int xb, int yb, int wb, int hb,
	int nx, int ix, float x1, float x2,
	int ny, int iy, float y1, float y2,
	int *nxb, int *ixb, float *x1b, float *x2b,
	int *nyb, int *iby, float *y1b, float *y2b);
static unsigned char *newInterpBytes (int n1in, int n2in, unsigned char *bin,
	int n1out, int n2out, int newInterpBytes);
void xMouseLoc(Display *dpy, Window win, XEvent event, int style, Bool show,
	int x, int y, int width, int height,
	float x1begb, float x1endb, float x2begb, float x2endb);
void xMousePrint(XEvent event, int style, FILE *mpicksfp,
	int x, int y, int width, int height,
	float x1begb, float x1endb, float x2begb, float x2endb);
/* JG... */
void intl2b_block(int nxin, float dxin, float fxin,
				  int nyin, float dyin, float fyin, unsigned char *zin,
				  int nxout, float dxout, float fxout,
				  int nyout, float dyout, float fyout, unsigned char *zout);
/* .... JG */



#ifdef TEST
main()
{
	int nxin,nyin,nxout,nyout,ixout,iyout;
	float dxin,fxin,dyin,fyin,fxout,dxout,fyout,dyout;
	unsigned char zin[2][2],zout[4][4];
	
	zin[0][0] = 41;		zin[0][1] = 99;
	zin[1][0] = 99;		zin[1][1] = 99;
	nxin=2;  dxin=1.0;  fxin=0.0;
	nyin=2;  dyin=1.0;  fyin=0.0;
	nxout=4;  dxout=dxin*(nxin-1)/(nxout-1);  fxout=0.0;
	nyout=4;  dyout=dyin*(nyin-1)/(nyout-1);  fyout=0.0;
	intl2b(nxin,dxin,fxin,nyin,dyin,fyin,&zin[0][0],
		nxout,dxout,fxout,nyout,dyout,fyout,&zout[0][0]);
	for (iyout=0; iyout<nyout; iyout++)
		for (ixout=0; ixout<nxout; ixout++)
			printf("zout[%d][%d] = %d\n",
				iyout,ixout,zout[iyout][ixout]);
}
#endif
void   requestdoc(int flag);
void   pagedoc(void);
char * cwp_strdup(char *str);
void *alloc1 (size_t n1, size_t size);
float *alloc1float(size_t n1) ;
int *alloc1int(size_t n1);
void free1 (void *p); 
void free1float(float *p); 
void free1int(int *p); 
void scaxis (float x1, float x2, int *nxnum, float *dxnum, float *fxnum);
void strchop(char *s, char *t);
void warn(char *fmt, ...);
static void intl2bx(int, int*, int*, int,
	unsigned char[][256], unsigned char*, unsigned char*);
static void intl2by(int, int, int, unsigned char[][256],
	unsigned char*, unsigned char*, unsigned char*);

void intl2b(int nxin, float dxin, float fxin,
	int nyin, float dyin, float fyin, unsigned char *zin,
	int nxout, float dxout, float fxout,
	int nyout, float dyout, float fyout, unsigned char *zout);

int xargc;  char **xargv;

int
main (int argc,char **argv)
{
	int n1,n2,n1tic,n2tic,
		i1,i2,grid1,grid2,style,
		n1c,n2c,i1beg,i1end,i2beg,i2end,i1c,i2c,
		nz,iz,i1step,i2step,verbose,
		xbox,ybox,wbox,hbox,
		xb,yb,wb,hb,
		x,y,width,height,
		i,j,nx,ny,nxb,nyb,ixb,iyb,
		imageOutOfDate,winwidth=-1,winheight=-1,
		showloc=0,
		legend,lwidth,lheight,lx,ly; /* BEREND */

	int blockinterp=0; /* JG */

	unsigned long nfloats;
	float labelsize,titlesize,perc,clip,bperc,wperc,bclip,wclip,
		d1,f1,d2,f2,*z,*temp,zscale,zoffset,zi,
		x1beg,x1end,x2beg,x2end,
		x1min,x1max,x2min,x2max,
		d1num,f1num,d2num,f2num,
		x1begb,x1endb,x2begb,x2endb,blank;
	unsigned char *cz,*czp,*czb,*czbp,*czbi=NULL;
	char *label1="",*label2="",*title="",*windowtitle="ximage",
		*units="", *legendfont="times_roman10",
		*labelfont="Erg14",*titlefont="Rom22",
		*styles="seismic",*grid1s="none",*grid2s="none",
		*labelcolor="blue",*titlecolor="red",
		*gridcolor="blue",*cmap="",keybuf[256],*mpicks;
	FILE *infp=stdin, *mpicksfp;
	Display *dpy;
	Window win;
	XEvent event;
	KeySym keysym;
	XComposeStatus keystat;
	XImage *image=NULL;
	XImage *image_legend=NULL; /* BEREND */
	unsigned char *data_legend; /* BEREND */
	GC gci;
	int scr;
	unsigned long black,white,pmin,pmax;

        float **x1curve,**x2curve;
        int curve,*npair,ncurvecolor;
        char **curvefile,**curvecolor=NULL;
        FILE *curvefp;

	/* initialize getpar */
	initargs(argc,argv);
	requestdoc(1);

	/* get parameters describing 1st dimension sampling */
	if (!getparint("n1",&n1))
		warn("Must specify number of samples in 1st dimension!\n");
	d1 = 1.0;  getparfloat("d1",&d1);
	f1 = 0.0;  getparfloat("f1",&f1);
	x1min = (d1>0.0)?f1:f1+(n1-1)*d1;
	x1max = (d1<0.0)?f1:f1+(n1-1)*d1;

	/* get parameters describing 2nd dimension sampling */
	if (!getparint("n2",&n2)) {
		if (fseeko(infp, (off_t) 0,SEEK_END)!=0)
			warn("must specify n2 if in a pipe!\n");
		nfloats = (int) (ftello(infp)/( (off_t) sizeof(float)));
		fseeko(infp, (off_t) 0,SEEK_SET);
		n2 = (int) (nfloats/n1);
	}
	d2 = 1.0;  getparfloat("d2",&d2);
	f2 = 0.0;  getparfloat("f2",&f2);
	x2min = (d2>0.0)?f2:f2+(n2-1)*d2;
	x2max = (d2<0.0)?f2:f2+(n2-1)*d2;

	/* set up file to save mouse picks */
	if (!getparstring("mpicks", &mpicks)) mpicks = "/dev/tty";
	mpicksfp = fopen(mpicks, "w");

	/* set up curve plotting */
	if ((curve=countparval("curve"))!=0) {
		curvefile=(char**)alloc1(curve,sizeof(void*));
		getparstringarray("curve",curvefile);
		if ((x1curve=(float**)malloc(curve*sizeof(void*)))==NULL)
			warn("Could not allocate x1curve pointers\n");
		if ((x2curve=(float**)malloc(curve*sizeof(void*)))==NULL)
			warn("Could not allocate x2curve pointers\n");
                npair=alloc1int(curve);
                getparint("npair",npair);
	} else {
		npair=(int *)NULL;
		curvefile=(char **)NULL;
		x1curve=(float **)NULL;
		x2curve=(float **)NULL;
	}
	if ((ncurvecolor=countparval("curvecolor"))<curve) {
		curvecolor=(char**)alloc1(curve,sizeof(void*));
		if (!getparstringarray("curvecolor",curvecolor)) {
			curvecolor[0]=(char *)cwp_strdup("blue\0");
			ncurvecolor=1;
		}
		for (i=ncurvecolor; i<curve; i++)
			curvecolor[i]=(char *)cwp_strdup(curvecolor[ncurvecolor-1]);
	} else if( ncurvecolor ) {
		curvecolor=(char**)alloc1(ncurvecolor,sizeof(void*));
		getparstringarray("curvecolor",curvecolor);
	}
	for (j=0; j<curve; j++) {
		curvefp=fopen(curvefile[j],"r");
		x1curve[j]=alloc1float(npair[j]);
		x2curve[j]=alloc1float(npair[j]);
		for (i=0; i<npair[j]; i++) {
			fscanf(curvefp,"%f",&x1curve[j][i]);
			fscanf(curvefp,"%f",&x2curve[j][i]);
		}
		fclose(curvefp);
	}

	/* read binary data to be plotted */
	nz = n1*n2;
	z = alloc1float(nz);
	if (fread(z,sizeof(float),nz,infp)!=nz)
		warn("error reading input file");

	/* if necessary, determine clips from percentiles */
	if (getparfloat("clip",&clip)) {
		bclip = clip;
		wclip = -clip;
	}
	if ((!getparfloat("bclip",&bclip) || !getparfloat("wclip",&wclip)) &&
		!getparfloat("clip",&clip)) {
		perc = 100.0;  getparfloat("perc",&perc);
		temp = alloc1float(nz);
		for (iz=0; iz<nz; iz++)
			temp[iz] = z[iz];
		if (!getparfloat("bclip",&bclip)) {
			bperc = perc;	getparfloat("bperc",&bperc);
			iz = (nz*bperc/100.0);
			if (iz<0) iz = 0;
			if (iz>nz-1) iz = nz-1;
			qkfind(iz,nz,temp);
			bclip = temp[iz];
		}
		if (!getparfloat("wclip",&wclip)) {
			wperc = 100.0-perc;  getparfloat("wperc",&wperc);
			iz = (nz*wperc/100.0);
			if (iz<0) iz = 0;
			if (iz>nz-1) iz = nz-1;
			qkfind(iz,nz,temp);
			wclip = temp[iz];
		}
		free1float(temp);
	}
	verbose = 1;  getparint("verbose",&verbose);
	if (verbose) warn("bclip=%g wclip=%g",bclip,wclip);

	/* get colormap specification */
	if (!(getparstring("cmap",&cmap))) {
		cmap = (char *)alloc1(5,1);
		sprintf(cmap,"%s","gray");
	}
	
	/* get interpolation style JG */
	if (getparint("blockinterp", &blockinterp)) blockinterp=1;

	/* get legend specs BEREND */
	legend = 0; getparint("legend", &legend); /* BEREND */
	getparstring("units", &units); /* BEREND */
	getparstring("legendfont", &legendfont);     /* BEREND */

	blank = 0; getparfloat("blank",&blank);

	/* get axes parameters */
	xbox = 50; getparint("xbox",&xbox);
	ybox = 50; getparint("ybox",&ybox);
	wbox = 550; getparint("wbox",&wbox);
	hbox = 700; getparint("hbox",&hbox);

	/* legend dimensions */
	if (!(getparint("lwidth",&lwidth)))	lwidth = 16;
	if (!(getparint("lheight",&lheight)))	lheight = hbox/3;
	if (!(getparint("lx",&lx)))	lx = 3;
	if (!(getparint("ly",&ly)))	ly = (hbox-lheight)/3;

	x1beg = x1min; getparfloat("x1beg",&x1beg);
	x1end = x1max; getparfloat("x1end",&x1end);
	d1num = 0.0; getparfloat("d1num",&d1num);
	f1num = x1min; getparfloat("f1num",&f1num);
	n1tic = 1; getparint("n1tic",&n1tic);
	getparstring("grid1",&grid1s);
	if (STREQ("dot",grid1s)) grid1 = DOT;
	else if (STREQ("dash",grid1s)) grid1 = DASH;
	else if (STREQ("solid",grid1s)) grid1 = SOLID;
	else grid1 = NONE;
	getparstring("label1",&label1);
	x2beg = x2min; getparfloat("x2beg",&x2beg);
	x2end = x2max; getparfloat("x2end",&x2end);
	d2num = 0.0; getparfloat("d2num",&d2num);
	f2num = 0.0; getparfloat("f2num",&f2num);
	n2tic = 1; getparint("n2tic",&n2tic);
	getparstring("grid2",&grid2s);
	if (STREQ("dot",grid2s)) grid2 = DOT;
	else if (STREQ("dash",grid2s)) grid2 = DASH;
	else if (STREQ("solid",grid2s)) grid2 = SOLID;
	else grid2 = NONE;
	getparstring("label2",&label2);
	getparstring("labelfont",&labelfont);
	labelsize = 18.0; getparfloat("labelsize",&labelsize);
	getparstring("title",&title);
	getparstring("titlefont",&titlefont);
	titlesize = 24.0; getparfloat("titlesize",&titlesize);
	getparstring("style",&styles);
	if (STREQ("normal",styles)) style = NORMAL;
	else style = SEISMIC;
	getparstring("titlecolor",&titlecolor);
	getparstring("labelcolor",&labelcolor);
	getparstring("gridcolor",&gridcolor);
	getparstring("windowtitle",&windowtitle);

	/* adjust x1beg and x1end to fall on sampled values */
	i1beg = NINT((x1beg-f1)/d1);
	i1beg = MAX(0,MIN(n1-1,i1beg));
	x1beg = f1+i1beg*d1;
	i1end = NINT((x1end-f1)/d1);
	i1end = MAX(0,MIN(n1-1,i1end));
	x1end = f1+i1end*d1;

	/* adjust x2beg and x2end to fall on sampled values */
	i2beg = NINT((x2beg-f2)/d2);
	i2beg = MAX(0,MIN(n2-1,i2beg));
	x2beg = f2+i2beg*d2;
	i2end = NINT((x2end-f2)/d2);
	i2end = MAX(0,MIN(n2-1,i2end));
	x2end = f2+i2end*d2;