dualplot.cpp

//**************************************************************************//
//* FRIDA: fast reliable interactive data analysis                         *//
//* dualplot.cpp: different mechanisms for screen and paper output         *//
//* (C) Joachim Wuttke 1990-, v2(C++) 2001-                                *//
//* http://frida.sourceforge.net                                           *//
//**************************************************************************//

#include <stdlib.h>
#include <iostream>
#include <math.h>

#include "mystd.h"
#include "range.h"
#include "grid.h"
#include "plotaux.h"
#include "coord.h"
#include "readln.h"
#include "asi.h"
#include "dualplot.h"

// local set-up:
#define GNUP    // Gnuplot
// #define PLUT // GNU plotutils (too low level; buffering problem)
// #define PLPL // PlPlot        (only Tek works, in X11 clipping problem)
#define WUPS    // my own postscript driver

// local implementation:
#ifdef GNUP
#include <sys/types.h>
#include <sys/stat.h>
#include <unistd.h>
#include <fcntl.h>
#endif

#ifdef PLUT
#include <plotter.h>
#endif

#ifdef PLPL
#include <plplot/plConfig.h>
#include <plplot/plplotP.h>
#endif

#ifdef WUPS
#define WUPS_GBUF "/tmp/rda.ps"
#endif

#define PLOTDTYP double
#define PSMAX 10



namespace NPlot {
	CAxis X, Y;

	uint maxpoints = 12000;

#ifdef PLUT
	Plotter *plut;
#endif
#ifdef GNUP
	int gp_fifo;
	void gp_write(string);
	int gp_fno;
	string gp_fnames;
#endif
#ifdef WUPS
	string ps_header;
	string ps_tmp;
	FILE *psout;
	uint  ps_fnum; // file
	uint  ps_snum; // scan
	uint  ps_pnum; // scan with pstyle
	uint  ps_cnum; // scan with cstyle
	void ps_ticktack(uint ntack, double *tack, int ntpt, double *ticklim,
			 CAxis *A);
	vector<string> ps_Doc;
#endif
};


void CAxis::Ask(const string& quest)
{
	string def, in, in1, in2;
	CRange newR;
	while(1) {
		def = (log ? "log " : "") + R.str();
		// out = quest + " [" + def + "] ?"; so gehts leider(?) nicht 
		in = sask(quest.c_str(), def);

		if        (in==def)
			return;

		mystd::string_extract_word(in, &in1, &in2);
		if        (in1=="l") {
			SetLog(!log);
			in = in2;
		} else if (in1=="i") {
			SetLog(0);
			in = in2;
		} else if (in1=="g") {
			SetLog(1);
			in = in2;
		} 

		if (in!="") {
			newR = CRange(in);
			if (!newR.valid)
				continue;
			if (log && newR.inf<=0) {
				printf("! log scale requires range above 0\n");
				continue;
			}
			R = newR;
			return;
		}
	}
}

void CAxis::SetLog(const uint val)
{
	if (log!=val) {
		log = val;
		CRange RX = AlternateR;
		AlternateR = R;
		R  = RX;
	}
}

#ifdef WUPS
double CAxis::pc(double v) // value -> plot-coordinate
{
	return PSMAX * R.relval(log, v);
}
#endif

void NPlot::Open(void)
{
#ifdef GNUP
	// Start gnuplot. Use input redirection so that gnuplot receives
        // commands from a gp_fifo which is created here.

    string fn_fifo = string("/tmp/gnuplot-") + getenv( "LOGNAME" );

	system( ("rm -f "+fn_fifo).c_str() );
	if (mkfifo( fn_fifo.c_str(), 0666 )) {
		printf("SEVERE SYSTEM ERROR cannot make fifo %s"
			" - will not be able to print\n", fn_fifo.c_str() );
		return;
	}

	system( ("gnuplot -noraise < " + fn_fifo + " &").c_str() );

	// we use open instead of fopen or ofstream,
        // because we need non-blocking mode.
       	if (!(gp_fifo = open(fn_fifo.c_str(), O_WRONLY))) {
		printf("SEVERE SYSTEM ERROR cannot open gp_fifo"
			" - will not be able to print\n");
		return;
	}
	fcntl(gp_fifo,F_SETFL,O_NONBLOCK);

	// gp_write(string("plot 1 2 3")); // test
	gp_write(string("set terminal x11"));
#endif

#ifdef PLUT
	// die folgenden Initialisierungen scheinen keinen Effekt zu haben:
	Plotter::parampl("VANISH_ON_DELETE", (void*)"no");
	Plotter::parampl("DISPLAY", (void*)":0.0");
	Plotter::parampl("USE_DOUBLE_BUFFERING", (void*)"no");

	plut = new XPlotter(cin, cout, cerr);

	if (plut->openpl() < 0) {
		printf("SEVERE SYSTEM ERROR "
		       "cannot open plotutils::Plotter::plut"
			" - will not be able to print\n");
		return;
	}
	plut->fspace(-.08, -.08, 1.04, 1.04);
#endif

#ifdef PLPL
	plsdev("xwin");
//	plsdev("xterm");
//	plsdev("tk");
	plscol0(0,255,255,255);
	plscol0(1,0,0,0);
	plSetOpt("np", ""); // no pause between plots
        //	plSetOpt("geometry", "600x400+340+5"); 
//	plSetOpt("debug", ""); 
	cout << "DEBUG P1\n";
	plinit();
	cout << "DEBUG P2\n";
	// plfont(2);
	// pltext(); // back to Tektronix text window
	cout << "DEBUG P3\n";
#endif

#ifdef WUPS
	ps_tmp = string("/tmp/rda-") + getenv("LOGNAME") + ".ps";
	ps_fnum=0;
#endif
}
		
void NPlot::Close(void)
{
#ifdef PLPL
	cout << "DEBUG PEND\n";
	plend();
#endif
}

void NPlot::Clear(void)
{
	// cout << "DEBUG dp Clear beg\n";
#ifdef GNUP
	gp_fno = 0;
	gp_fnames = string("");
#endif

#ifdef PLUT
#endif

#ifdef PLPL
	// new pl window:
	// plgra(); // go to Tektronix graph window
	pladv(0);
	plvpor(0.15, 0.85, 0.1, 0.9);
#endif

#ifdef WUPS
	// new wups file:
	if (NRead::wofmac("\\pshead", &ps_header)) {
		printf("! please define \\pshead\n");
		return;
	}
	string cmd;
	cmd = "cp -f " + ps_header + " " + ps_tmp;
	if (system(cmd.c_str())) {
		printf("! command [%s] failed\n", cmd.c_str());
		return;
	} else if (!(psout = fopen(ps_tmp.c_str(), "r"))) {
		printf("copy has not created %s\n", ps_tmp.c_str());
		return;
	}
	fclose(psout);
	if (!(psout = fopen(ps_tmp.c_str(), "a"))) {
		printf("cannot write (append) to %s\n", ps_tmp.c_str());
		return;
	}
	fprintf(psout, "\n%% output created by rda\n");
	ps_snum = 0;
	ps_pnum = 0;
	ps_cnum = 0;
	ps_Doc.clear();
#endif
}

void NPlot::Multiask(void)
{
  	maxpoints = iask("Max # points in plot", maxpoints);
	X.force = bask("Force x points into frame", X.force);
	Y.force = bask("Force y points into frame", Y.force);
}

void NPlot::Box()
{
#ifdef GNUP
	string whichlog="", cmd;
	if (X.log) whichlog += "x";
	if (Y.log) whichlog += "y";
	if (whichlog=="")
		cmd = "set nologscale";
	else
		cmd = "set logscale " + whichlog;
	gp_write(cmd);
#endif

#ifdef PLUT
	plut->flinewidth(0.001);
	plut->pencolorname("black");
	plut->fbox(0., 0., 1., 1.);
	plut->fbox(0., 0., 1., 1.);
	plut->flushpl();
#endif

#ifdef PLPL
	plwind(X.R.inf, X.R.sup, Y.R.inf, Y.R.sup);
	plcol0(1);
	plbox("bcnst", 0.0, 0, "bcnstv", 0.0, 0);
	// plbox("bcfghlnst", 0.0, 0, "bcghnstv", 0.0, 0);

	// draw labels:
    /*
	plcol0(1);
	// plmtex("b", 3.2, 0.5, 0.5, "Frequency");
	// plmtex("t", 2.0, 0.5, 0.5, "Single Pole Low-Pass Filter");
	plcol0(2);
	// plmtex("l", 5.0, 0.5, 0.5, "Amplitude (dB)");
     */
#endif

#ifdef WUPS
	// wups:
	fprintf(psout, "\n%d %f %f xSetCoord\n", X.log, X.R.inf, X.R.sup);
	fprintf(psout, "%d %f %f ySetCoord\n", Y.log, Y.R.inf, Y.R.sup);
	fprintf(psout, "%% %d %f %f %d zSetCoord\n\n", 0, 0., 0., 0);

	int ntack, ntpt;
	double *tack, ticklim[2];

	fprintf(psout, "\n/xPlotFrame {");
  	if(plotaux::calc_ticks(X.log, X.R.inf, X.R.sup, 
			       &ntack, &tack, &ntpt, ticklim)) {
  		printf("PROGRAM ERROR/ Ticks\n");
		return;
	}
	ps_ticktack(ntack, tack, ntpt, ticklim, &X);
	free(tack);
	fprintf(psout, "   {(%s)}\n", X.C.ps_str().c_str());
	fprintf(psout, 
  "   schoolstyle 1 eq {\n"
  "      0 10   0  0 wy  0  90 ArrAxx Axx xTacC %% xNumL %% low x axis\n"
  "   } {\n"
  "      0 10   0  0     0  90 OneAxx Axx Tic xTacL xNumL %% low x axis\n"
  "      0 10   0 10     0 270 OneAxx Axx Tic xTacH       %% top x axis\n"
  "      xCL\n"
  "   } ifelse\n"
  "} def\n");

	fprintf(psout, "\n/yPlotFrame {");
  	if(plotaux::calc_ticks(Y.log, Y.R.inf, Y.R.sup, 
			       &ntack, &tack, &ntpt, ticklim)) {
  		printf("PROGRAM ERROR/ Ticks\n");
		return;
	}
	ps_ticktack(ntack, tack, ntpt, ticklim, &Y);
	free(tack);
	fprintf(psout, "   {(%s)}\n", Y.C.ps_str().c_str());
	fprintf(psout, 
  "   schoolstyle 1 eq {\n"
  "      0 10   0 wx 0  90   0 ArrAxx Axx yTacC %% yNumL %% left y axis\n"
  "   } {\n"
  "      0 10   0  0    90   0 OneAxx Axx Tic yTacL yNumL %% left y axis\n"
  "      0 10  10  0    90 180 OneAxx Axx Tic yTacH %% yNumH %% right y axis\n"
  "      yCL\n"
  "   } ifelse\n"
  "} def\n");

	fprintf(psout, "\n%% modeDD\nplotbefore\n");
#endif
}

void NPlot::Line(const int lstyle, 
		 const CGrid x, const vector<double> y,
		 const vector<double>* z,
		 const string xco, const string yco)
{
	PLOTDTYP *xp, *yp;
	uint np=0, nxl=0, nxh=0, nyl=0, nyh=0;
	uint n=x.size();
	if (n!=y.size()) {
		printf("PROG ERR NPLot::Line x.size<>y.size\n");
		return;
	}
	xp = (PLOTDTYP *) malloc(sizeof(PLOTDTYP)*n);
	yp = (PLOTDTYP *) malloc(sizeof(PLOTDTYP)*n);
	for (uint i=0; i<n; i++) {
		if(!X.force && !X.R.contained(x[i])) {
			if (x[i]<=X.R.inf) nxl++;
			if (x[i]>=X.R.sup) nxh++;
			continue;
		}
		if(!Y.force && !Y.R.contained(y[i])) {
			if (y[i]<=Y.R.inf) nyl++;
			if (y[i]>=Y.R.sup) nyh++;
			continue;
		}
		if(np==maxpoints && np<n) {
			Tekx();
			printf("reached maxpoints at %g\n", x[i]);
		}
		xp[np] = x[i];
		yp[np] = y[i];
		np++;
		if(np>maxpoints) i += maxpoints;
	}
	if (np==0) {
		Tekx();
		printf("no points in data range:\n"
		      "of %d points are\n %d left and %d right of the x range,"
		      "\n %d below and %d above of the y range\n", 
		       n, nxl, nxh, nyl, nyh);
		return;
	}

#ifdef GNUP
	// save to file:
	char gp_fnam[40];
	sprintf( gp_fnam, "/tmp/%s-%03d.gnu", getenv("LOGNAME"), gp_fno++);
	if (gp_fnames!="") gp_fnames += ", ";
	gp_fnames += string("\"") + gp_fnam + "\"" + 
//		string(" title \"") + xco + string (" -> ") 
//		                    + yco + string ("\"");
		" notitle";
	if (lstyle<0)
		gp_fnames += " with lines";
	FILE *gp_fd;
	if (!(gp_fd = fopen(gp_fnam, "w"))) {
		printf( "cannot save gnuplot data to %s", gp_fnam );
		return;
	} 
	for (uint i=0; i<np; i++)
		fprintf(gp_fd, "%20.12g %20.12g\n", xp[i], yp[i]);
	fclose(gp_fd);
	
	// plot command:
	char aux[80];
	sprintf(aux, "[%20.12g:%20.12g] [%20.12g:%20.12g] ", 
		X.R.inf, X.R.sup, Y.R.inf, Y.R.sup);
	gp_write(string("plot ") + aux + gp_fnames);
#endif

#ifdef PLUT
#endif

#ifdef PLPL
	plcol0(9);
	plpoin(np, xp, yp, 1);
#endif

#ifdef WUPS
	if (!psout) {
		printf("internal file psout currently not open\n");
		return;
	}
	fprintf(psout, "\n%3u [", ++ps_snum);
	for (uint i=0; i<z->size(); i++) 
		fprintf(psout, " %12g", (*z)[i]);
	fprintf(psout, " ] zValues\n");

	if (lstyle>=0)
		fprintf(psout, "%2d pstyle", ++ps_pnum);
	else
		fprintf(psout, "%2d cstyle", 1);
	fprintf(psout, " %% (%s -> %s)\n", xco.c_str(), yco.c_str());
	for (uint i=0; i<np; i++) {
		fprintf(psout, "%7.3f%7.3f%7.3f t",
			X.pc(xp[i]), Y.pc(yp[i]), 0.0);
		if      (i==0)    fprintf(psout, "i");
		else if (i==np-1) fprintf(psout, "f");
		else              fprintf(psout, " ");
		fprintf(psout, "%% %- 12.6g wx %- 12.6g wy\n", xp[i], yp[i]);
	}
#endif

	free(xp);
	free(yp);
}

void NPlot::Doc (vector<string> lDoc)
{
#ifdef WUPS
	for (uint i=0; i<lDoc.size(); i++)
		ps_Doc.push_back(lDoc[i]);
#endif	
}

void NPlot::Save(void)
{
#ifdef WUPS
        if (!psout) {
	        printf("internal file psout currently not open\n");
		return;
	}
	fclose(psout);

	string fpref = NRead::exemac("\\psdir");
	if (fpref!=string("")) fpref += string("/");

	FILE *pssav;
	char outf[20];
	string flong, cmd;
	while(1) {
		// construct not-yet-used file name:
		if (ps_fnum>=999) {
			printf ("graph file number overflow\n");
			return;
		}
		sprintf(outf, "%sl%d.%s", fpref.c_str(), ++ps_fnum, "psa");
		if (!(pssav = mystd::glob_fopen(outf, "", "", "r")))
			break; // legal exit
		fclose(pssav);
	}
	printf("save plot in %s\n", outf);
	cmd = string("cp ") + ps_tmp + " " + outf + "\n";
	if (system(cmd.c_str())) {
		printf ("cannot copy to %s\n", outf);
		return;
	}
	if (!(pssav = mystd::glob_fopen(outf, "", "", "r", &flong))) {
		printf("copy has not created file %s\n", outf);
		return;
	}
	fclose(pssav);
	if (!(pssav = fopen(flong.c_str(), "a+"))) {
		printf ("cannot write (append) to  file %s\n", flong.c_str());
		return;
	}

	fprintf(pssav, "\nblack 0 -4 13 newlist\n");
	for (uint i=0; i<ps_Doc.size(); i++) 
		fprintf(pssav, "{(%s)} infline\n", ps_Doc[i].c_str());
	
	fprintf(pssav, 
		"\n{(%s)}  /filename exch def 10 -2.8 18 showfilename\n\n"
		"EndFrame\n", outf);

	fclose(pssav);

	if (!(psout = fopen(ps_tmp.c_str(), "a"))) {
		printf("cannot reopen internal file psout [%s]\n", ps_tmp.c_str());
		return;
	}
#else
	cout << "graphic save desactivated\n";
#endif
#ifdef PLPL
//	plend(); // TEST
#endif
}

void NPlot::Tekx(void)
{
#ifdef PLPL
	// pltext(); // back to Tektronix text window
#endif
}

void NPlot::Dialog(void)
{
#ifdef GNUP
	string cmd;
	cout << "entering mygnuplot - to leave, type q\n";
	while(1) {
		cmd = sask("mygnuplot> ");
		if (cmd.substr(0,1)=="q") return;
		gp_write(cmd);
	}
#else
	cout << "Graphic dialog is implemented only for Gnuplot\n";
#endif
}
#ifdef GNUP
void NPlot::gp_write(string in)
{
	// cout << "DEBUG gp_write " << in << "\n";
	string out = in + "\n";
	write(gp_fifo, out.c_str(), out.size());
}
#endif

#ifdef WUPS
void NPlot::ps_ticktack(uint ntack, double *tack, int ntpt, double *ticklim,
			CAxis *A)
{
	uint i;
	fprintf(psout, "[\n");
	if (A->log && ( tack[0]<1e-3 || tack[ntack-1]>1e3 )) {
	    for (i=0; i<ntack; i++) {
			fprintf(psout, "   %9.6f {(10)(%d)sp()} %%{(%g)}\n", 
				A->pc(tack[i]), mystd::nint(log10(tack[i])),
				(float) tack[i]);
	    }
	} else {
		for (i=0; i<ntack; i++) {
			fprintf(psout, "   %9.6f {(%g)}\n", 
				A->pc(tack[i]), (float) tack[i]);
		}
	}
	fprintf(psout, "   ] SetTacVec\n");
	fprintf(psout, "   %g %g %d %d SetTicVec%s\n", 
		A->pc(ticklim[0]), A->pc(ticklim[1]), ntack+2, ntpt,
		(A->log? "Log" : "Lin"));
}
#endif