dualplot.cpp

//**************************************************************************************************
//*  FRIDA: fast reliable interactive data analysis
//*  (C) Joachim Wuttke 1990-, v2(C++) 2001-
//*  http://apps.jcns.fz-juelich.de/frida
//**************************************************************************************************

//! \file  dualplot.cpp
//! \brief Collection NPlot of plot frames CPlot.

#include <iostream>
#include <string>
#include <vector>
#include <cstdio>
#include <cmath>
#include <cstring>
#include <fcntl.h>
#include <errno.h>
#include <sys/stat.h>
#include <boost/format.hpp>

#include "../trivia/file_ops.hpp"
#include "../trivia/string_convs.hpp"

#include "dualplot.hpp"

using std::string;
using std::vector;
using std::cout;
using boost::format;

#define S(a) triv::strg((a))

static const int mLin = 80; // max num of chars in PostScript line


//! Constructor for plot window: setup for gnuplot and postscript.

CPlot::CPlot( int _iPlot, bool _logx, bool _logy ) :
    iPlot( _iPlot ), X( "x", _logx ), Y( "y", _logy ),
    maxpoints(20000), with_errors(true), equipoints(49), refine(true)
{
    // == Initialization for Gnuplot ==

    // Create a named pipe (FIFO) that will transmit our commands to Gnuplot.
    string fn_fifo = string("/tmp/gnuplot-") + getenv( "LOGNAME" );
    triv::system( "rm -f "+fn_fifo );
    if (mkfifo( fn_fifo.c_str(), 0666 ))
        throw "SYSTEM ERROR cannot make fifo "+fn_fifo+": will not be able to print";

    // Check that Gnuplot supports X11.
    string out = triv::system_read( "gnuplot -e 'help X11' 2>&1 < /dev/null" );
    if ( out.substr(0,5)=="Sorry" || out.length()<80 ) {
        cout << "Gnuplot seems not to support X11\n";
        exit(1);
    }

    // Start a Gnuplot that will read from our pipe.
    triv::system( "gnuplot -title "+S(iPlot)+" -noraise < "+fn_fifo+"&" );

    // Open our pipe so that we can write commands to it
    //   (we use 'open' instead of 'fopen' or 'ofstream', because we need non-blocking mode).
    if ( ( gp_fifo = open( fn_fifo.c_str(), O_WRONLY ) )==-1 )
        throw "SYSTEM ERROR: cannot open fifo "+fn_fifo+" ("+strerror(errno)+")";
    fcntl( gp_fifo, F_SETFL, O_NONBLOCK );

    // Now the initialization _within_ Gnuplot.
    gp_write( string("set terminal x11") );

    // == Initialization for PostScript ==
    ps_fnum=0;
}


//! Clear plot frame and buffer.

void CPlot::clear_frame()
{
    // clear gnuplot tmp files
    gp_fno = 0;
    gp_fnames = "";

    // reset buffers for postscript output
    string cmd;
    ps_accu.clear();
    ps_accu.push_back( "\n%% output created by frida2\n");
    ps_snum = 0;
    ps_pnum = 0;
    ps_cnum = 0;
    ps_Doc.clear();
}


//! Change one setup parameter per command.

void CPlot::set_aux( const string& cmd )
{
    if ( cmd=="gxl" ) {
        X.set_log( !X.logflag );
        printf( "set x %s\n", X.logflag ? "log" : "lin" ); }
    else if ( cmd=="gxl+" )
        X.set_log( true );
    else if ( cmd=="gxl-" )
        X.set_log( false );
    else if ( cmd=="gyl" ) {
        Y.set_log( !Y.logflag );
        printf( "set y %s\n", Y.logflag ? "log" : "lin" ); }
    else if ( cmd=="gyl+" )
        Y.set_log( true );
    else if ( cmd=="gyl-" )
        Y.set_log( false );
    else if ( cmd=="gxf" ) {
        X.force = !X.force;
        printf( "force x %s\n", X.force ? "on" : "off" ); }
    else if ( cmd=="gxf+" )
        X.force = true;
    else if ( cmd=="gxf-" )
        X.force = false;
    else if ( cmd=="gyf" ) {
        Y.force = !Y.force;
        printf( "force y %s\n", Y.force ? "on" : "off" ); }
    else if ( cmd=="gyf+" )
        Y.force = true;
    else if ( cmd=="gyf-" )
        Y.force = false;
    else if ( cmd=="ge" )
        with_errors = !with_errors;
    else if ( cmd=="ge+" )
        with_errors = true;
    else if ( cmd=="ge-" )
        with_errors = false;
    else
        throw "unknown command " + cmd;
}


//! Plot coordinate frame (axes, ticks, labels).

void CPlot::plot_frame( const string& xlabel, const string& ylabel )
{
    gp_write( "set nologscale" );
    string whichlog="";
    if (X.logflag) whichlog += "x";
    if (Y.logflag) whichlog += "y";
    if (whichlog!="")
        gp_write( "set logscale " + whichlog );

    // wups:
    char outlin[ mLin ];
    snprintf( outlin, mLin, "\n%i %g %g xSetCoord\n",
              X.logflag, X.inf, X.sup );
    ps_accu.push_back( outlin );
    snprintf( outlin, mLin, "%i %g %g ySetCoord\n",
              Y.logflag, Y.inf, Y.sup );
    ps_accu.push_back( outlin );
    snprintf( outlin, mLin, "%% %i %g %g %i zSetCoord\n\n",
              0, 0., 0., 0 );
    ps_accu.push_back( outlin );

    int ntpt;
    double ticklim[2];
    vector<double> Tacks;
    ps_accu.push_back( "\n/xPlotFrame {\n" );
    if ( X.logflag && X.inf<= 0 )
        throw "BUG: x log incompatible with limits " + X.str();
    X.calc_ticks( Tacks, &ntpt, ticklim );
    ps_ticktack( Tacks, ntpt, ticklim, &X);
    snprintf( outlin, mLin-4, "  {(%s", xlabel.c_str() );
    strncat( outlin, ")}\n", mLin );
    ps_accu.push_back( outlin );
    ps_accu.push_back( "   0 10   0  0     0  90 "
                       "OneAxx Axx Tic Tac xNumL %% low x axis\n" );
    ps_accu.push_back( "   0 10   0 10     0 270 "
                       "OneAxx Axx Tic Tac       %% top x axis\n" );
    ps_accu.push_back( "  xCL\n" );
    ps_accu.push_back( "} def\n" );

    ps_accu.push_back( "\n/yPlotFrame {\n" );
    if ( Y.logflag && Y.inf<= 0 )
        throw "BUG: y log incompatible with limits " + Y.str();
    Y.calc_ticks( Tacks, &ntpt, ticklim );
    ps_ticktack( Tacks, ntpt, ticklim, &Y);
    snprintf( outlin, mLin-4, "   {(%s", ylabel.c_str() );
    strncat( outlin, ")}\n", mLin );
    ps_accu.push_back( outlin );
    ps_accu.push_back( "   0 10   0  0    90   0 "
                       "OneAxx Axx Tic Tac yNumL %% left y axis\n" );
    ps_accu.push_back( "   0 10  10  0    90 180 "
                       "OneAxx Axx Tic Tac       %% right y axis\n" );
    ps_accu.push_back( "  yCL\n" );
    ps_accu.push_back( "} def\n" );
    ps_accu.push_back( "\n%% modeDD\nplotbefore\n" );
}


//! Plot one spectrum.

void CPlot::add_spec( bool as_line, bool new_style, int style_no,
                     const vector<double>& xp,
                     const vector<double>& yp, const vector<double>& dyp,
                     const vector<string>& z,
                     const string& xco,
                     const string& yco,
                     const string& info
    )
{
    static const int mColor = 6;
    static int color[mColor] = { 0x880000, 0x008800, 0x000088,
                                 0x006666, 0x660066, 0x666600 };
    // Checks:
    int np=xp.size();
    if ( !np )
        throw S("invalid call to CPLot::addSpec: no data points");
    if ( np!=yp.size() )
        throw S("invalid call to CPLot::addSpec: x.size<>y.size");

    // Prepare for live display, to be shown by showSpecs():
    string gp_fnam = str( format( "/tmp/%s-%i-%03i.gnu" )
                          % getenv("LOGNAME") % iPlot % gp_fno++ );
    if (gp_fnames!="") gp_fnames += ", ";
    gp_fnames += string("\"") + gp_fnam + "\" notitle";
    if ( as_line )
        gp_fnames += str( format( " with lines lt 1 lc rgb \"#%6x\"" )
                          % color[ style_no % mColor ] );
    else if( with_errors && dyp.size() )
        gp_fnames += " with errorbars";
    FILE *gp_fd;
    if (!(gp_fd = fopen(gp_fnam.c_str(), "w")))
        throw "cannot save gnuplot data to " + gp_fnam;
    int nout = 0;
    try {
        for (int i=0; i<np; i++){
            if( std::isinf(xp[i]) || std::isinf(yp[i]) )
                throw "Data point number " + S(i) + " is invalid: x=" +
                    S(xp[i]) + ", y=" + S(yp[i]);
            if( xp[i]<X.inf || xp[i]>X.sup )
                throw "CPlot::addSpec: x["+S(i)+"]="+S(xp[i])+" out of range";
            if( yp[i]<Y.inf || yp[i]>Y.sup )
                throw "CPlot::addSpec: y["+S(i)+"]="+S(yp[i])+" out of range";
            if( with_errors && dyp.size() )
                fprintf(gp_fd, "%16.8g %16.8g %16.8g\n", xp[i], yp[i], dyp[i] );
            else
                fprintf(gp_fd, "%16.8g %16.8g\n", xp[i], yp[i]);
            nout++;
        }
    } catch ( string &s ) {
        fclose(gp_fd);
        throw s;
    } catch ( ... ) {
        fclose(gp_fd);
        throw "BUG: unexpected exception type";
    }
    fclose(gp_fd);
    if( !nout )
        throw "no points in frame: " + info;

    // Postscript copy:
    char outlin[ mLin ];
    if ( new_style ) {
        snprintf( outlin, mLin, "\n%3u [", ++ps_snum );
        ps_accu.push_back( outlin );
        for (int i=0; i<z.size(); i++){
            snprintf( outlin, mLin, " %s", z[i].c_str() );
            ps_accu.push_back( outlin );
        }
        snprintf( outlin, mLin, " ] zValues\n" );
        ps_accu.push_back( outlin );
        if ( as_line )
            snprintf( outlin, mLin, "%2i cstyle", style_no );
        else
            snprintf( outlin, mLin, "%2i pstyle", style_no );
        ps_accu.push_back( outlin );
        snprintf( outlin, mLin-2, " %% (%s -> %s)", xco.c_str(), yco.c_str() );
        strncat( outlin, "\n", mLin );
        ps_accu.push_back( outlin );
    } else {
        ps_accu.push_back( "\n" );
    }
    for (int i=0; i<np; i++) {
        snprintf( outlin, mLin,
                  "%6.3f %6.3f %6.3f t%c %% %13.7g wx %13.7g wy\n",
                  X.pc(xp[i]), Y.pc(yp[i]),
                  dyp.size() ? Y.pcerr(yp[i],dyp[i]) : 0,
                  i==0 ? 'i' : i==np-1 ? 'f' : ' ',
                  xp[i], yp[i] );
        ps_accu.push_back( outlin );
    }
}


//! Live display as prepared by addSpec.

void CPlot::show_specs()
{
    if ( gp_fnames!="" )
        gp_write( "plot " +
                  str( format( "[%12.8g:%12.8g] [%12.8g:%12.8g] " )
                       % X.inf % X.sup % Y.inf % Y.sup ) +
                  gp_fnames );
}


//! Add documentation line to postscript output.

void CPlot::doc_TxLine( const string& line )
{
    ps_Doc.push_back( "  {("+line+")} TxLine" );
}


//! Add documentation line explaining a plot symbol to postscript output.

void CPlot::doc_PtTxLine( const string& line, int num )
{
    ps_Doc.push_back( "  " + S(num) + " {("+line+")} PtTxLine" );
}


//! Add documentation line explaining a curve style to postscript output.

void CPlot::doc_CvTxLine( const string& line, int num )
{
    ps_Doc.push_back( "  " + S(num) + " {("+line+")} CvTxLine" );
}


//! Write buffered plot to postscript file.

void CPlot::write_postscript(
    const string& ps_outdir, const string& ps_head, const string& ps_dict )
{
    // construct output file name:
    FILE *pssav;
    string cmd, outf;
    while(1) {
        if (ps_fnum>=999)
            throw S("graph file number overflow");
        outf = triv::wordexp_unique(
            ps_outdir + str( format( "l%i" ) % ++ps_fnum ) + "." + ( ps_dict=="" ? "psa" : "ps" ) );
        if( !triv::file_exists( outf.c_str() ) )
            break; // legal exit
    }
    cout << "save plot in " << outf << "\n";

    // copy headers to output file:
    cmd = string("cat ") +
        ps_dict + " " + // ps_dict may be ""
        ps_head + " > " +
        outf + "\n";
    triv::system( cmd );

    // append specific output to output file:
    if ( !(pssav = fopen( outf.c_str(), "a+" )) )
        throw "cannot append contents to file " + outf;
    for( string lin: ps_accu ){
        // fprintf does not work here because output line may contain "%"
        fwrite( lin.c_str(), 1, lin.size(), pssav );
    }

    // additional output (do not append this to ps_accu to allow
    // further incrementation of ps_accu):
    fprintf( pssav, "\n{ black 0 -4 13 1.65 NewList\n" );
    for ( string lin: ps_Doc )
        fprintf( pssav, "%s\n", lin.c_str() );
    fprintf( pssav, "} oooinfo 1 eq { exec } { pop } ifelse\n" );

    fprintf( pssav,
             "\n{(%s)}  /filename exch def 10 -2.8 18 showfilename\n\n"
            " EndFrame\n", outf.c_str() );

    // output completed:
    fclose(pssav);
}


//! Info line to characterize this plot window.

string CPlot::info() const
{
    string ret;
    ret  =   "x: " + X.info();
    ret += "  y: " + Y.info();
    return ret;
}


//! Send one line to gnuplot fifo.

void CPlot::gp_write( const string& in )
{
    string out = in + "\n";
    // cout << "monitor gnuplot driver: '" << out << "'\n";
    if( write( gp_fifo, out.c_str(), out.size() ) <= 0 )
        throw S("could not write to gp_fifo");
}


//! Format ticks and tacks for postscript file.

void CPlot::ps_ticktack( const vector<double>& Tacks, int ntpt,
                         const double *ticklim, const CAxis *A )
{
    char outlin[ mLin ];
    int i, ntack;
    ntack = Tacks.size();
    if (ntack > 0 ) {
        ps_accu.push_back( "  [\n" );
        if (A->logflag && ( Tacks[0]<1e-3 || Tacks[ntack-1]>1e3 )) {
            for (i=0; i<ntack; i++) {
                snprintf( outlin, mLin,
                          "   %9.6f {(10)(%i)sp()} %%{(%g)}\n",
                          A->pc(Tacks[i]), (int)(log10(Tacks[i])), (float) Tacks[i]);
                ps_accu.push_back( outlin );
            }
        } else {
            for (i=0; i<ntack; i++) {
                snprintf( outlin, mLin, "   %9.6f {(%g)}\n",
                          A->pc(Tacks[i]), (float) Tacks[i]);
                ps_accu.push_back( outlin );
            }
        }
        ps_accu.push_back( "  ] SetTacVec\n" );
    }
    snprintf( outlin, mLin, "  %g %g %i %i SetTicVec%s\n",
              A->pc(ticklim[0]), A->pc(ticklim[1]), ntack+2, ntpt,
              (A->logflag? "Log" : "Lin"));
    ps_accu.push_back( outlin );
}