pscube.c

su 的源代码库

    bclip = clip;
    wclip = -clip;
  }
  if ((!getparfloat("bclip",&bclip) || !getparfloat("wclip",&wclip)) &&
      !getparfloat("clip",&clip)) {
    perc = 100.0;  getparfloat("perc",&perc);
    temp = ealloc1float(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);
  }
  if (!getparint("verbose",&verbose)) verbose = 1;
  if (verbose) warn("bclip=%g wclip=%g",bclip,wclip);
  
  /* get scaled sampling intervals */
  if (!getparfloat("d1s",&d1s)) d1s = 1.0;
  if (!getparfloat("d2s",&d2s)) d2s = 1.0;
  if (!getparfloat("d3s",&d3s)) d3s = 1.0;
  d1s = fabs(d1s);  d1s *= d1;
  d2s = fabs(d2s);  d2s *= d2;
  d3s = fabs(d3s);  d3s *= d3;
  
  /* get projection angle, convert to radians */
  if(!getparfloat("angle",&angle)) angle = 45.0;
  angle = MAX(angle,0.00001);
  angle = MIN(angle,90.0);
  angle *= PI/180.0;
  
  /* get axes parameters */
  if(!getparfloat("size1",&size1)) size1 = 4.0;
  if(!getparfloat("size2",&size2)) size2 = 4.0;
  if(!getparfloat("size3",&size3)) size3 = 3.0;
  if (!getparfloat("xbox",&xbox)) xbox = 1.5;
  if (!getparfloat("ybox",&ybox)) ybox = 1.5;
 
  /* compute extreme values */
  x1min = (d1>0.0)?f1:f1+(n1-1)*d1;
  x1max = (d1<0.0)?f1:f1+(n1-1)*d1;
  x2min = (d2>0.0)?f2:f2+(n2-1)*d2;
  x2max = (d2<0.0)?f2:f2+(n2-1)*d2;
  x3min = (d3>0.0)?f3:f3+(n3-1)*d3;
  x3max = (d3<0.0)?f3:f3+(n3-1)*d3;
  
  /* get axis1 parameters */
  x1beg = x1min;
  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);
  
  /* get axis2 parameters */
  x2beg = x2min; getparfloat("x2beg",&x2beg);
  x2end = x2max; 
  d2num = 0.0; getparfloat("d2num",&d2num);
  f2num = x2min; 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);
  
  /* get axis3 parameters */
  x3beg = x3min; 
  x3end = x3max; getparfloat("x3end",&x3end);
  d3num = 0.0; getparfloat("d3num",&d3num);
  f3num = x3min; getparfloat("f3num",&f3num);
  n3tic = 1; getparint("n3tic",&n3tic);
  getparstring("grid3",&grid3s);
  if (STREQ("dot",grid3s)) grid3 = DOT;
  else if (STREQ("dash",grid3s)) grid3 = DASH;
  else if (STREQ("solid",grid3s)) grid3 = SOLID;
  else grid3 = NONE;
  getparstring("label3",&label3);
  
  /* get additional font parameters */
  getparstring("labelfont",&labelfont);
  labelsize = 18.0; getparfloat("labelsize",&labelsize);
  getparstring("title",&title);
  getparstring("titlefont",&titlefont);
  titlesize = 24.0; getparfloat("titlesize",&titlesize);
  getparstring("titlecolor",&titlecolor);
  getparstring("axescolor",&axescolor);
  getparstring("gridcolor",&gridcolor);
  style = SEISMIC;
  
 
  /* adjust x1beg and x1end to fall on sampled values */
  i1beg = NINT((x1beg-f1)/d1);
  i1beg = MAX(0,MIN(n1,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;
  
  /* adjust x3beg and x3end to fall on sampled values */
  i3beg = NINT((x3beg-f3)/d3);
  i3beg = MAX(0,MIN(n3-1,i3beg));
  x3beg = f3+i3beg*d3;
  i3end = NINT((x3end-f3)/d3);
  i3end = MAX(0,MIN(n3-1,i3end));
  x3end = f3+i3end*d3;

  /* allocate space for image bytes */
  n1c = 1+abs(i1end-i1beg);
  n2c = 1+abs(i2end-i2beg);
  n3c = 1+abs(i3end-i3beg);
  czfront = ealloc1(n1c*n2c,sizeof(char));
  czside = ealloc1(n1c*n3c,sizeof(char));
  cztop = ealloc1(n2c*n3c,sizeof(char));
  
  /* compute conversion constants */
  zscale = (wclip!=bclip)?255.0/(wclip-bclip):1.0e10;
  zoffset = -bclip*zscale;
  i1step = (i1end>i1beg)?1:-1;
  i2step = (i2end>i2beg)?1:-1;
  i3step = (i3end>i3beg)?1:-1;
  
  /* convert front data to be imaged into unsigned characters */
  czp = czfront;
  for (i2c=0,i2=i2beg; i2c<n2c; i2c++,i2+=i2step) {
    for (i1c=0,i1=i1beg; i1c<n1c; i1c++,i1+=i1step) {
      zi = zoffset+zfront[i1+i2*n1]*zscale;
      if (zi<0.0) zi = 0.0;
      if (zi>255.0) zi = 255.0;
      *czp++ = (unsigned char)zi;
    }
  }
  
  /* convert side data to be imaged into unsigned characters */
  czp = czside;
  for (i3c=0,i3=i3beg; i3c<n3c; i3c++,i3+=i3step) {
    for (i1c=0,i1=i1beg; i1c<n1c; i1c++,i1+=i1step) {
      zi = zoffset+zside[i1+i3*n1]*zscale;
      if (zi<0.0) zi = 0.0;
      if (zi>255.0) zi = 255.0;
      *czp++ = (unsigned char)zi;
    }
  }

  /* convert top data to be imaged into unsigned characters */
  czp = cztop;
  for (i3c=0,i3=i3beg; i3c<n3c; i3c++,i3+=i3step) {
    for (i2c=0,i2=i2beg; i2c<n2c; i2c++,i2+=i2step) {
      zi = zoffset+ztop[i2+i3*n2]*zscale;
      if (zi<0.0) zi = 0.0;
      if (zi>255.0) zi = 255.0;
      *czp++ = (unsigned char)zi;
    }
  }
  free1float(z);
  
  /* determine sampling after scaling */
  n1s = MAX(1,NINT(1+(n1c-1)*d1/d1s));
  d1s = (n1s>1)?d1*(n1c-1)/(n1s-1):d1;
  n2s = MAX(1,NINT(1+(n2c-1)*d2/d2s));
  d2s = (n2s>1)?d2*(n2c-1)/(n2s-1):d2;
  n3s = MAX(1,NINT(1+(n3c-1)*d3/d3s));
  d3s = (n3s>1)?d3*(n3c-1)/(n3s-1):d3;
  
  /* if necessary, interpolate front to scaled sampling intervals */
  if (n1s!=n1c || n2s!=n2c) {
    szfront = ealloc1(n1s*n2s,sizeof(char));
    intl2b(n1c,d1,0.0,n2c,d2,0.0,czfront,
	   n1s,d1s,0.0,n2s,d2s,0.0,szfront);
    free1(czfront);
  } else {
    szfront = czfront;
  }
  
  /* if necessary, interpolate side to scaled sampling intervals */
  if (n1s!=n1c || n3s!=n3c) {
    szside = ealloc1(n1s*n3s,sizeof(char));
    intl2b(n1c,d1,0.0,n3c,d3,0.0,czside,
	   n1s,d1s,0.0,n3s,d3s,0.0,szside);
    free1(czside);
  } else {
    szside = czside;
  }
  
  /* if necessary, interpolate top to scaled sampling intervals */
  if (n2s!=n2c || n3s!=n3c) {
    sztop = ealloc1(n2s*n3s,sizeof(char));
    intl2b(n2c,d2,0.0,n3c,d3,0.0,cztop,
	   n2s,d2s,0.0,n3s,d3s,0.0,sztop);
    free1(cztop);
  } else {
    sztop = cztop;
  }
  
  /* determine axes pads */
  p1beg = (x1end>x1beg)?-fabs(d1s)/2:fabs(d1s)/2;
  p1end = (x1end>x1beg)?fabs(d1s)/2:-fabs(d1s)/2;
  p2beg = (x2end>x2beg)?-fabs(d2s)/2:fabs(d2s)/2;
  p2end = (x2end>x2beg)?fabs(d2s)/2:-fabs(d2s)/2;
  p3beg = (x3end>x3beg)?-fabs(d3s)/2:fabs(d3s)/2;
  p3end = (x3end>x3beg)?fabs(d3s)/2:-fabs(d3s)/2;
  
  /* get legend specs BEREND, Schoenfelder */
  legend = 0; getparint("legend", &legend); /* BEREND, Schoenfelder */
  getparstring("units", &units); /* BEREND, Schoenfelder */
  getparstring("legendfont", &legendfont);     /* BEREND, Schoenfelder */
  
  /* Get or calc legend parameters */
  /* Legend min and max: Calc from data read in */
  if (legend) {
    for (lz=0;lz<nz;lz++) {
      lmin=FMIN(lmin,z[lz]);
      lmax=FMAX(lmax,z[lz]);
    }
    if (verbose==2) warn("lmin=%g lmax=%g",lmin,lmax);
  
    lbeg = lmin; if (getparfloat("lbeg",&lbeg)) lbegsup=1;
    lend = lmax; if (getparfloat("lend",&lend)) lendsup=1;
    /* Change wclip,bclip to be inside legend range */
    wclip = FMAX(lbeg,wclip); /* [wclip,bclip] has to be in [lbeg,lend] */
    bclip = FMIN(lend,bclip);
    if (lbegsup!=1) { /* Add white and black areas to show possible clipping */ 
      float rangeperc=(bclip-wclip)/20.;
      lbeg=wclip-rangeperc;
    }
    if (lendsup!=1) {
      float rangeperc=(bclip-wclip)/20.;
      lend=bclip+rangeperc;
    }
    
    lfnum = lmin; getparfloat("lfnum",&lfnum);
    
    getparstring("lstyle",&lstyles);
    if (STREQ("vertright",lstyles))
      lstyle = VERTRIGHT;
    else if (STREQ("horibottom",lstyles))
      lstyle = HORIBOTTOM;
    /* legend dimensions (BEREND), Schoenfelder */
    lwidth = 0.1 ;lheight = size1+sin(angle)*size3/2;
    if (lstyle==HORIBOTTOM) {
      lwidth=size2+cos(angle)*size3/1.2 ;lheight = 0.24;
    }
    getparfloat("lwidth",&lwidth);
    getparfloat("lheight",&lheight);
    
    lx=.8;ly = ybox+(size1+sin(angle)*size3-lheight)/2;
    if (lstyle==VERTRIGHT) {
      lx=xbox+size2+cos(angle)*size3+0.1;
    } else if (lstyle==HORIBOTTOM) {
      lx=xbox+(size2+cos(angle)*size3-lwidth)/2.0;ly = 1.0;
    }
    getparfloat("lx",&lx);
    getparfloat("ly",&ly);
    
    getparstring("lgrid",&lgrids);
    if (STREQ("dot",lgrids))
      ugrid = DOT;
    else if (STREQ("dash",lgrids))
      ugrid = DASH;
    else if (STREQ("solid",lgrids))
      ugrid = SOLID;
    else
      ugrid = NONE;
  }
  

  if (legend) {
    /* Make legend color values */
    int lll=0,lcount,perc5=13,ilbeg,ilend; /* color scale */
    if (lbegsup!=1) {
      ln+=perc5; /* white area */
    }
    if (lendsup!=1) {
      ln+=perc5; /* black area */
    }
    data_legend = ealloc1(ln,sizeof(char));
    if (lbegsup!=1) {
      for (lll=0;lll<perc5;lll++) data_legend[lll]=(char) 255; /* white area */
    }