Forked from
mlz / Frida
2348 commits behind the upstream repository.
-
Wuttke, Joachim authoredWuttke, Joachim authored
edif.cpp 38.12 KiB
//**************************************************************************//
//* FRIDA: fast reliable interactive data analysis *//
//* edif.cpp: edit on-line files (text_edit, dir, make, plot, ...) *//
//* (C) Joachim Wuttke 1990-, v2(C++) 2001- *//
//* http://frida.sourceforge.net *//
//**************************************************************************//
#include <stdlib.h>
#include <math.h>
#include <iostream>
#include <algorithm>
#include "mystd.h"
#include "olm.h"
#include "expr.h"
#include "asi.h"
#include "dualplot.h"
#include "edif.h"
#include "index.h"
using namespace std;
namespace NEdif {
CParam* RFind(CEle* e, CCoord OneR);
}
//! Display list of online files.
void NEdif::Dir(void)
{
CEle *e;
NOlm::IterateEle iter;
iter.All();
size_t npts;
string txtpts;
while((e=iter())) {
npts = e->nPts();
if (npts)
txtpts = strg(npts);
else
txtpts = "?";
char ztext[17];
if( e->P()->nZ() ){
if ( e->P()->ZCo.size()==1 )
sprintf( ztext, "%16s",
e->P()->ZCo[0].str().substr(0,16).c_str() );
else if( e->P()->ZCo.size()==2 )
sprintf( ztext, "%16s",
( e->P()->ZCo[0].str().substr(0,10) + ", " +
e->P()->ZCo[1].str().substr(0,4) ).c_str() );
else if( e->P()->ZCo.size()==3 )
sprintf( ztext, "%16s",
( e->P()->ZCo[0].str().substr(0,5) + ", " +
e->P()->ZCo[1].str().substr(0,4) + ", " +
e->P()->ZCo[2].str().substr(0,3) ).c_str() );
else
sprintf( ztext, "%16s",
( e->P()->ZCo[0].str().substr(0,4) + ", " +
e->P()->ZCo[1].str().substr(0,3) + ", " +
e->P()->ZCo[2].str().substr(0,2) + " &c"
).c_str() );
} else {
sprintf( ztext, "" );
}
printf("%3d%1s%-16s %5d %-16s %6s %-12s %-12s\n",
iter.SelNo(),
e->P()->as_on_disk ? "=" : " ",
e->P()->name.substr(0,16).c_str(), e->nScan(),
ztext,
txtpts.c_str(),
e->P()->xco.str().substr(0,12).c_str(),
e->P()->yco.str().substr(0,12).c_str()
);
}
}
//! Display list of z-values and x-ranges for scans in selected files.
void NEdif::DirZ(void)
{
NOlm::IterateEle iter;
const CEle *e;
const COld *d;
const COlc *c;
while((e=iter())) {
cout << "file " << iter.SelNo() << ":\n";
printf(" ");
c = e->C();
d = e->D();
for (size_t m=0; m<e->P()->ZCo.size(); m++) {
printf(" %12s", e->P()->ZCo[m].str().c_str());
}
printf( "%8s %38s\n", "#pts.", "x" );
for (size_t j=0; j<e->nScan(); j++) {
printf("%3zu", j);
for (size_t m=0; m<e->P()->ZCo.size(); m++) {
printf(" %12.6g", (*(e->Z(j)))[m]);
}
if (d)
printf( "%12.6g (%4zu) %12.6g\n",
d->VS[j].x[0], d->VS[j].size(),
d->VS[j].x[d->VS[j].size()-1] );
else
printf("\n");
}
}
}
//! Display z coordinates (names and units) for selected files.
void NEdif::ListZ(void)
{
NOlm::IterateEle eiter;
CEle *e;
size_t k, nz, i;
string act;
vector<CCoord> ZCo;
if (!(e=eiter())) return;
// Get current coordinates:
ZCo = e->P()->ZCo;
nz = ZCo.size();
if (nz)
for (i=0; i<nz; ++i)
printf(" %zu: %s\n", i, ZCo[i].str().c_str());
else
printf(" there are no z coordinates\n");
// Check other input files:
while((e=eiter())) {
k = eiter.SelNo();
if ( e->P()->ZCo.size()!=nz )
throw "file " + strg(k) + " has different number of z coordinates";
for (i=0; i<nz; ++i)
if ( e->P()->ZCo[i]!=ZCo[i] )
throw "file " + strg(k) +
" has different z coordinate " + strg(i);
}
}
//! Display numeric parameters for selected files.
void NEdif::DirNumpar(void)
{
NOlm::IterateEle iter;
const CEle *e;
while((e=iter())) {
cout << "file " << iter.SelNo() << ":\n";
printf(" ");
for (size_t m=0; m<e->P()->RPar.size(); m++) {
cout << " " << e->P()->RPar[m].str() << "\n" ;
}
}
}
//! Interactive modification of numeric parameters.
//! QUESTIONABLE. Why not use generic expressions and specific commands?
void NEdif::EditNumpar(void)
{
NOlm::IterateEle eiter;
CEle *e;
vector<CCoord> AllR;
CParam *rp;
vector<uint> RStat; // 1: different values 2: partially undefined
vector<double> RVal;
string choice;
IndexSet Choice;
IndexSetIterator ISet;
CCoord newco;
double val;
uint i, ii;
bool first;
make_table:
// Collate RPar entries:
AllR.clear();
eiter.Reset();
while((e=eiter())) {
for (i=0; i<e->P()->RPar.size(); ++i) {
for (ii=0; ii<AllR.size(); ++ii)
if (e->P()->RPar[i].Co == AllR[ii])
goto found;
AllR.push_back(e->P()->RPar[i].Co);
found:
;
}
}
RStat.reserve(AllR.size());
RVal.reserve(AllR.size());
for (ii=0; ii<AllR.size(); ++ii) {
RStat[ii] = 0;
first = 1;
eiter.Reset();
while((e=eiter())) {
if ((rp=RFind(e, AllR[ii]))) {
if (first) {
RVal[ii] = rp->val; first = 0;
} else if (RVal[ii] != rp->val) {
RStat[ii] |= 1;
}
} else {
RStat[ii] |= 2;
}
}
}
// Show table:
for (ii=0; ii<AllR.size(); ++ii) {
printf(" %3u %.40s", ii, AllR[ii].str().c_str());
if (RStat[ii] & 1)
printf("<different values> ");
else
printf("%18.10g ", RVal[ii]);
if (RStat[ii] & 2)
printf("<partially undefined>\n");
else
printf("\n");
}
user_choice:
choice = wask("Modify r-pars (h=help) ?", 3, "q");
if (choice == "h") {
printf("Enter either list of r-par numbers"
" or a:add\n");
goto user_choice;
} else if (choice == "a") {
newco.Ask("New parameter");
val = dask("Value");
eiter.Reset();
while((e=eiter()))
e->P()->RPar.push_back(CParam(newco, val));
goto make_table;
} else if (choice=="" || choice=="q") {
return;
} else if (Choice.parse(choice, 0, AllR.size())) {
cout << mystd::BEL;
goto user_choice;
}
// It's a list
// User action on list of parameters:
ISet.reset(Choice);
while (ISet(&ii)) {
if (RStat[ii]) {
// Table:
/*
eiter.Reset();
while((e=eiter())) {
*/
}
action_on_par:
choice = wask(("Action on r-par[" + strg(ii) +"] <val,d,q,<>>"
).c_str());
if (choice=="h") {
printf(" q quit\n"
" d delete entry\n"
" <val> set uniform value\n"
" f set value per file\n");
goto action_on_par;
} else if (choice=="q") {
return;
} else if (choice=="") {
} else if (choice=="d") {
/* eiter.Reset();
while((e=eiter())) {
for (i=0; i<e->P()->RPar.size(); ++i)
if (e->P()->RPar[i]
== AllR[ii]) {
e->P()->RPar.erase(
printf("missing ...\n");
*/
} else if (sscanf(choice.c_str(), "%lg", &val)==1) {
eiter.Reset();
while((e=eiter())) {
if ((rp=RFind(e, AllR[ii])))
rp->val = val;
else
e->P()->RPar.push_back(
CParam(AllR[ii], val));
}
} else {
cout << mystd::BEL;
goto action_on_par;
} // action on rpar
} // list of rpars
goto make_table;
}
//! Auxiliary routine for EditNumpar.
CParam* NEdif::RFind(CEle* e, CCoord OneR)
{
for (size_t i=0; i<e->P()->RPar.size(); ++i)
if (e->P()->RPar[i].Co == OneR) return &(e->P()->RPar[i]);
return 0;
}
//! Display coordinate names and units for selected files.
void NEdif::DirCoord(void)
{
NOlm::IterateEle eiter;
CEle *e;
string out;
char aux[4];
while((e=eiter())) {
sprintf( aux, "%3d", eiter.SelNo() );
out = string( aux );
out += " x: " + e->P()->xco.str();
out += " y: " + e->P()->yco.str();
for (size_t i=0; i<e->P()->nZ(); i++) {
out += " z" + strg(i) + ": " + e->P()->ZCo[i].str() + "\n";
}
cout << out;
}
}
//! Modify file name
void NEdif::EditFNam()
{
NOlm::IterateEle eiter;
const CEle *ein;
COlo *f;
string fnam;
while( (ein=eiter()) ) {
f = ein->P();
fnam = wask( ("Rename " + f->name).c_str(), f->name );
if ( fnam != f->name ){ f->lDoc.push_back( "ef " + fnam + " # < " + f->name );
f->name = fnam;
}
}
}
//! Modify coordinate names and units.
void NEdif::EditCoord( string which )
{
NOlm::IterateEle eiter;
const CEle *ein; CEle eout;
COlo *fin; COlo *fout;
CRef cref( which );
if ( !cref.defined() ) {
printf( "! invalid coordinate %s\n", which.c_str() );
return;
}
CCoord new_co, old_co;
new_co.Ask( (which+" coordinate").c_str() );
// no-overwrite mode is broken since 28aug08
while( (ein=eiter()) ) {
fin = ein->P();
fout = fin;
if ( cref.var==CRef::_X ) {
old_co = fin->xco;
fout->xco = new_co;
} else if ( cref.var==CRef::_Y ) {
old_co = fin->yco;
fout->yco = new_co;
} else if ( cref.var==CRef::_Z ) {
if( cref.num >= fin->nZ() ){
printf( "! not all files in selection have %s, use oz+\n",
which.c_str() );
return;
}
old_co = fin->ZCo[cref.num];
fout->ZCo[cref.num] = new_co;
}
fout->lDoc.push_back( "ec"+string(which)+" "+new_co.str()+
" # old: " + old_co.str() );
}
}
//! Display documentation lines for selected files.
//! MISSING: Interactive modification of documentation.
void NEdif::EditDoc(void)
{
NOlm::IterateEle iter;
CEle *e;
while( (e=iter()) ){
COlo* f = e->P();
COlc* fc = e->C();
if( iter.index()>0 )
cout << "\n";
cout << iter.SelNo() << "\n";
for ( size_t i=0; i<f->lDoc.size(); i++ )
cout << f->lDoc[i] << "\n";
if( fc )
cout << fc->pInfoCat();
}
}
//! Save online file as external file (the work is done in COlo::Save).
void NEdif::Save( string fmt ){
COlo *f;
string fnam;
FILE *fd;
CEle *e;
NOlm::IterateEle iter;
while((e=iter())) {
if( !(f = e->P()) )
throw string( "SEVERE / f ?" );
string quest = "Save as (."+fmt+")";
fnam = wask( quest.c_str(), f->name);
if ( fnam==string("") ) return;
if ( fnam != f->name && fmt=="y08" ){
f->lDoc.push_back( "fs " + fnam + " # < " + f->name );
f->name = fnam;
}
fnam += "."+fmt;
if ( (fd = fopen(fnam.c_str(), "r")) ) {
fclose(fd);
if (!bask("Overwrite ?")) return;
}
if ( !(fd = fopen(fnam.c_str(), "w")) )
throw "! cannot write to file " + fnam;
e->Save( fd, fmt );
fclose(fd);
}
}
//! Create a new online file from interactive input.
void NEdif::ReadIn(void)
{
COld olf;
// *** set file name and coordinates ***
static string fnam;
olf.name = wask("Save as", fnam);
static CCoord xco, yco, zco;
olf.xco = xco.Ask("x", xco);
olf.yco = yco.Ask("y", yco);
cout << "if there is only one scan, answer \"\\\" to the following:\n";
zco.Ask("z", zco);
if(zco.defined()) olf.ZCo.push_back(zco);
// *** set input mode ***
static int mInpX = -1;
cout << "input modes:\n"
" (1) enter x pointwise\n"
" (2) set regular x values\n"
" (3) enter x-y pairs\n"
" (4) select x,y from multicolumn input\n"
" (5) select x,y,z from multicolumn input\n";
mInpX = iask( "Option", mInpX );
if( mInpX==5 && !zco.defined() )
throw string( "need z coordinate for this option" );
static int imcx, imcy, imcz;
int mInpY=-1; // ask
switch(mInpX) {
case 0:
return; case 1: case 2:
break;
case 3:
imcx = 0;
imcy = 1;
mInpY = 0;
break;
case 4: case 5:
imcx = iask("x from column #", imcx);
if (imcx<0) return;
imcy = iask("y from column #", imcy);
if (imcy<0) return;
mInpY = 0;
break;
default:
throw string( "invalid option" );
}
if (mInpY<0) {
cout << "input modes:\n"
" (1) enter y pointwise\n"
" (2) set regular y values\n"
" (3) select y from multicolumn input\n";
static int mInpY_in = -1;
mInpY_in = iask( "Option", mInpY_in );
if ( mInpY_in<1 || mInpY_in>3 )
throw string( "invalid option" );
mInpY = mInpY_in;
}
switch(mInpY) {
case 0: case 1: case 2:
break;
case 3:
imcy = iask("y from column #", imcy);
if (imcy<0) return;
break;
default:
throw string( "PROGRAM ERROR in Make: invalid switch(mInpY)" );
}
static int mInpZ;
if ( mInpX==5 ){
mInpZ = 3;
} else if (zco.defined()) {
cout << "at present, we allow only _one_ z value per scan\n";
cout << "input modes:\n"
" (1) enter z per spectrum\n"
" (2) select z from multi-z header line\n";
mInpZ = iask( "Option", mInpZ );
if ( mInpZ<0 || mInpZ>2 )
throw string( "invalid option" );
} else {
mInpZ = 0;
}
switch (mInpZ) {
case 0: case 1:
break;
case 2:
if (! (mInpX==4 || mInpY==3) )
throw string( "invalid combination of options:"
" y is not read from multicolumn format" );
break;
case 3:
imcz = iask( "z0 from column", imcz );
if( imcz<0 ) return;
break;
default:
throw string( "PROGRAM ERROR in Make: invalid switch(mInpZ)" ); }
static int imcy_delta=0;
static uint imcy_last;
cout << "DEBUG "<< mInpZ << "\n";
if (mInpZ==2) {
imcy_delta = iask("step in number of y column", imcy_delta);
if (imcy_delta)
imcy_last = iask("last y column", imcy_last);
} else {
imcy_delta = 0;
}
// *** loop over scans ***
CScan S;
uint n = 0;
uint imcm;
string quest, line;
double val, zold;
vector<double> linv;
int ns = 0;
while(1) {
S.Clear();
switch(mInpZ) {
case 1:
line = sask("Enter z value", "");
// cout << "DEBUG: line [" << line << "]\n";
if (mystd::str2vec(line, &linv, 2) < 0) {
cout << "-> invalid line [" << line << "]\n";
goto store;
}
// printf ("DEBUG: size = %d\n", linv.size());
if (linv.size()<1) {
cout << "found no z in line [" << line << "] => assuming end-of-input\n";
goto store;
}
// printf ("DEBUG: val = %g\n", linv[0]);
val = linv[0];
break;
case 2:
line = sask("Enter header line with z values", "");
if (mystd::str2vec(line, &linv, imcy+1)) {
cout << "-> bad input\n";
goto store;
}
if (linv.size()<=imcy) {
cout << "-> not enough entries in header\n";
goto store;
}
val = linv[imcy];
break;
}
if (mInpZ<=2)
S.z.push_back(val);
switch(mInpX) {
case 1:
cout << "Enter x values (or empty line to exit):\n";
while (1) {
quest = "[" + strg(n) + "] ";
line = sask(quest.c_str(), "");
if (sscanf(line.c_str(), "%lg", &val)!=1) break;
S.x.push_back(val);
n++;
}
break; case 2:
//S.x.Ask("x");
//n = S.x.size();
throw string( "case 2 broken" );
break;
case 3: case 4: case 5:
cout<<"Enter data row (or empty line to exit):\n";
imcm = (imcx>imcy) ? imcx : imcy;
if( mInpY==5 && imcz>imcm )
imcm = imcz;
while (1) {
quest = "[" + strg(n) + "] ";
line = sask(quest.c_str(), "");
if (line==string("")) break;
if (mystd::str2vec(line, &linv, imcm+1)) {
printf(
" could not extract value no. %u"
" from string [%s]\n", imcm, line.c_str());
continue;
}
if (linv.size()<=imcm) {
cout << "-> not enough entries in line\n";
continue;
}
if( mInpZ==3 ){
if( n==0 ){
S.z.push_back(linv[imcz]);
zold = linv[imcz];
} else if( linv[imcz]!=zold ){
olf.VS.push_back(S);
ns++;
S.Clear();
S.z.push_back(linv[imcz]);
zold = linv[imcz];
}
}
S.push_back(linv[imcx], linv[imcy]);
n++;
}
break;
}
if (n<=0) return;
switch(mInpY) {
case 0:
break;
case 1:
cout << "Enter " << n << " y values:\n";
S.y.resize(n);
for (uint i=0; i<n; i++) {
quest = "x[" + strg(i) + "] = " +
strg(S.x[i]) + " -> ";
do {
line = sask(quest.c_str());
} while (sscanf(line.c_str(), "%lg", &val)!=1);
S.y[i] = val;
}
break;
case 2:
throw string( "case 2 broken, ygrid undefined" );
//ygrid.Ask("y", n);
//ygrid.AsVec(&(S.y));
break;
case 3:
cout << "Enter " << n << " lines with y values:\n";
S.y.resize(n);
for (uint i=0; i<n; i++) {
quest = "x[" + strg(i) + "] =" +
strg(S.x[i]) + " -> ";
while (1) { line = sask(quest.c_str());
if (mystd::str2vec(line,
&linv, imcy+1)) {
printf(
" could not extract y value no. %u"
" from string [%s]\n", imcy, line.c_str());
continue;
}
if (linv.size()<=imcy) {
cout << "-> not enough entries"
" in line\n";
continue;
}
break;
}
S.y[i] = linv[imcy];
}
break;
}
olf.VS.push_back(S);
ns++;
if (!mInpZ || mInpZ==3) goto store;
if (imcy_delta) {
imcy += imcy_delta;
if (imcy_delta>0 && imcy>imcy_last ||
imcy_delta<0 && imcy<imcy_last) goto store;
}
}
// *** store new on-line file ***
store:
// cout << "DEBUG: store\n";
if (ns<=0) return;
// cout << "DEBUG: sel\n";
NOlm::SelNew();
// cout << "DEBUG: add\n";
NOlm::OloAdd(&olf);
// cout << "DEBUG: main\n";
}
//! Create a new online file from table.
void NEdif::ReadTab( string qualif )
{
// ** parse qualif **
if( qualif.find_first_not_of("hvcsm")!=string::npos )
throw string( "ReadTab: invalid qualifier" );
char dir = qualif[0];
if( !( dir=='h' || dir=='v' || dir=='c' ) )
throw string( "ReadTab: missing qualifier h or v" );
bool horizontal = dir=='h';
bool choosecol = dir=='c';
static int iycol=0;
bool fromscript = qualif.find('s')!=string::npos;
bool multiblock = qualif.find('m')!=string::npos;
// ** input from script or from file **
static string script;
vector<string> inFiles;
if( fromscript ){
inFiles.push_back( "/ram/tab" );
cout << "Script must write to "<<inFiles[0]<<"\n";
script = wask("Script", script);
mystd::system( script );
} else {
string fnames = wask( "Read tab from file(s)" );
// PROVIS, sask sobald glob besser
if( mystd::glob_file_list( fnames, "", "", &inFiles )<=0 ) throw "files " + fnames + " not found";
}
if( choosecol ){
printf( "at present, x:=i\n" );
iycol = iask( "Read y from column", iycol );
if( iycol<0 ) return;
}
FILE *fd;
string fdir, fshort, fext;
NOlm::SelNew();
for( uint iF=0; iF<inFiles.size(); ++iF ) {
if( !(fd = fopen(inFiles[iF].c_str(), "r")) )
throw string( "cannot open file " ) + inFiles[iF];
cout << ".. reading from " << inFiles[iF] << "\n";
COld olf;
if( fromscript ) {
olf.lDoc.push_back( "ft"+qualif+" " + script );
olf.name = script;
} else {
mystd::fname_divide( inFiles[iF], &fdir, &fshort, &fext);
olf.lDoc.push_back( "ft"+qualif+" " + inFiles[iF] );
olf.name = fshort;
if( choosecol ){
olf.lDoc.push_back( "y from column " + strg(iycol) );
olf.name += "_" + strg(iycol);
}
}
// ** file-level settings **
if( choosecol ){
olf.xco = CCoord("i", "");
olf.yco = CCoord("y"+strg(iycol), "");
} else {
olf.xco = CCoord("x", "");
olf.yco = CCoord("y", "");
}
// *** read input ***
string lin, s1, s2;
vector<double> dat, zdat;
int n_in = -1; // line number
int nblock = 0;
int nline = -1; // line number within block (starting out of block)
int nz = 0;
int nzdat = 0;
double val;
CScan S;
while( mystd::freadln(fd, &lin)>-1 ) {
++n_in;
if( lin.substr(0,4)=="#rpa" ){
mystd::string_extract_word( lin.substr(4), &s1, &s2 );
if( mystd::any2dbl( s2, &val ) )
printf( "invalid parameter line [%s]\n", lin.c_str() );
else
olf.RPar.push_back( CParam( s2, "", val) );
continue;
}
if( lin[0]=='#' )
continue;
if( nline==-1 ) { // start of block
if( !horizontal && nblock!=0 && S.size()>0 )
olf.VS.push_back(S);
if ( multiblock ) {
if( mystd::str2vec(lin, &zdat, 0, false) )
throw "invalid header line [" + lin +
"] (or legitimate break ?)";
if( nblock==0 ){ // first header nzdat = zdat.size();
nz += nzdat;
for( uint iz=0; iz<nzdat; ++iz )
olf.ZCo.push_back( CCoord("z"+strg(iz), "") );
} else if( zdat.size() != nzdat )
throw "line " + strg(n_in) +
": header has length " + strg(zdat.size()) +
" instead of " + strg(nzdat);
}
if ( nblock==0 && horizontal ) {
nz += 1;
olf.ZCo.push_back( CCoord("line", "") );
}
if( !horizontal ){
S.Clear();
S.z.resize( nz );
for( uint iz=0; iz<nzdat; ++iz )
S.z[iz] = zdat[iz];
}
nline = 0;
++nblock;
if( multiblock ) // this was a header line
continue;
}
if( lin=="" && nline!=-1 ) { // end-of-block
nline = -1;
++nblock;
continue;
}
// regular data line
if( mystd::str2vec(lin, &dat, 0, false) )
throw "cannot parse line "+strg(nline)+" in block "+
strg(nblock)+" ["+lin+"]";
if(dat.size()==0)
throw "no data in line "+strg(nline)+" in block "+
strg(nblock)+" ["+lin+"]";
if( horizontal ) { // new scan for every line
S.Clear();
S.z.resize( nz );
for( uint iz=0; iz<nzdat; ++iz )
S.z[iz] = zdat[iz];
S.z[nz-1] = nline;
for( int i=0; i<dat.size(); ++i )
S.push_back( (double)i, dat[i] );
olf.VS.push_back(S);
} else { // vertical
if( choosecol ){
if( iycol>=dat.size() )
throw "line "+strg(n_in)+" (line "+strg(nline)+
" in block "+strg(nblock)+") ["+lin+
"] contains only "+strg(dat.size())+" values; "+
"cannot read y from col "+strg(iycol);
S.push_back( (double)nline, dat[iycol] );
} else {
if( dat.size()!=2 )
throw "line "+strg(n_in)+" (line "+strg(nline)+
" in block "+strg(nblock)+") ["+lin+
"] contains "+strg(dat.size())+" values; "+
"at present, exactly 2 are required";
S.push_back( dat[0], dat[1] );
}
}
++nline;
} // end of file input loop
if( !horizontal && S.size()>0 ) olf.VS.push_back(S);
if( !(olf.nScan()) )
throw string( "no input lines" );
if( olf.nScan()==1 ){
olf.ZCo.clear();
olf.VS[0].z.clear();
}
NOlm::OloAdd(&olf);
}
}
//! Create a new online file containing a regular x,z grid; let y=0.
void NEdif::MakeGrid(void)
{
COld olf;
static int ni=1, nj=1;
static string fnam;
int niIn = iask("Number of points per scan", ni);
if (niIn<1) return;
ni = niIn;
int njIn = iask("Number of scans", nj);
if (njIn<1) return;
nj = njIn;
fnam = "grid"+strg(ni);
if (nj>1) fnam += "x"+strg(nj);
olf.name = wask("Save as", fnam);
// *** set coordinates ***
olf.xco = CCoord("x", "");
olf.yco = CCoord("y", "");
if (nj>1)
olf.ZCo.push_back(CCoord("z", ""));
CScan S;
// *** loop over scans ***
for (int j=0; j<nj; ++j) {
S.Clear();
S.x.resize(ni);
for( uint i=0; i<ni; ++i )
S.x[i] = ((double)i)/(ni+1);
S.y.clear();
S.y.resize( ni, 0. );
if (nj>1) S.z.push_back( (double)j );
olf.VS.push_back(S);
}
NOlm::SelNew();
NOlm::OloAdd(&olf);
}
//! Plot scan.
void NEdif::Plot( CPlot *plot, bool add )
{
static CList JSel;
CEle *e; size_t k, j;
const COld *fd;
const COlc *fc;
CScan S;
static int pstyle = 0, cstyle = 0;
static string jSel = "";
NOlm::JSelAsk( "Plot which scans", &jSel, &JSel );
if (!add) {
pstyle = 0;
cstyle = 0;
// determine x-range ?
double inf, sup;
if (!plot->X.finite()) {
NOlm::IterateD fiter;
const vector<double> *vx;
bool first = true;
while((fd=fiter())) {
JSel.evaluate( 0, fd->nScan()-1 );
for( size_t iv=0; iv<JSel.size(); ++iv ) {
vx = &( fd->VS[JSel.V[iv]].x );
for( size_t i=0; i<vx->size(); ++i ){
if( plot->X.logflag && (*vx)[i]<=0 )
continue;
if( first ){
inf = (*vx)[i];
sup = (*vx)[i];
first = false;
} else {
if( (*vx)[i]<inf ) inf = (*vx)[i];
if( (*vx)[i]>sup ) sup = (*vx)[i];
}
}
}
}
if( first )
throw string( "x range is empty" );
plot->X.setLimits( inf, sup );
plot->X.round();
}
// determine y-range ?
if (!plot->Y.finite()) {
bool first=true;
NOlm::IterateEle eiter;
int npts = 987;
const CScan *s;
while( (e=eiter()) ){
k = eiter.SelNo();
fc = e->C();
fd = e->D();
JSel.evaluate( 0, e->nScan()-1 );
for( size_t iv=0; iv<JSel.size(); ++iv ){
size_t j = JSel.V[iv];
if( fd ){
s = &(fd->VS[j]);
} else {
S.x.resize( npts );
for ( size_t i=0; i<npts; i++ )
S.x[i] = plot->X.plotcoord2value( i/(npts-1.0) );
fc->T->tree_curve_val_vec( &(S.y), S.x, k, j );
s = &S;
}
for (size_t i=0; i<s->size(); i++) {
if( !plot->X.contained( s->x[i] ) )
continue;
if( plot->Y.logflag && s->y[i]<=0 ) continue;
if (first) {
inf=s->y[i], sup=s->y[i];
first = false;
} else {
if (s->y[i]<inf) inf=s->y[i];
if (s->y[i]>sup) sup=s->y[i];
}
}
}
}
if( first )
throw string( "y range is empty" );
plot->Y.setLimits( inf, sup );
plot->Y.round();
}
// build label: (zur Zeit nur vom ersten File; definiere +=)
plot->X.C = NOlm::MOM[*(NOlm::FSel.V.begin())].P()->xco;
plot->Y.C = NOlm::MOM[*(NOlm::FSel.V.begin())].P()->yco;
// draw new frame:
plot->clearFrame();
plot->plotFrame();
}
// plot data:
NOlm::IterateEle eiter;
while( (e=eiter()) ){
k = eiter.SelNo();
fc = e->C();
fd = e->D();
JSel.evaluate( 0, e->nScan()-1 );
for ( size_t iv=0; iv<JSel.size(); ++iv){
j = JSel.V[iv];
// cout << "DEBUG edif plot k="<<k<<", j="<<j<<"\n";
if ( fd ) {
const CScan *s;
s = &( fd->VS[j] );
size_t n=s->x.size();
if ( n!=s->y.size() )
throw string( "BUG: plot: x.size<>y.size" );
// Filter: x,y -> xp,yp if inside frame
size_t np=0, nxl=0, nxh=0, nyl=0, nyh=0;
vector<double> xp, yp;
for ( size_t i=0; i<n; i++ ) {
double x = s->x[i];
if( !plot->X.contained(x) ) {
if ( x<=plot->X.inf ){
x = plot->X.inf;
nxl++;
}
if ( x>=plot->X.sup ){
x = plot->X.sup;
nxh++;
}
if ( !plot->X.force )
continue;
}
double y = s->y[i];
if( !plot->Y.contained(y) ) {
if ( y<=plot->Y.inf ){
y = plot->Y.inf;
nyl++;
}
if ( y>=plot->Y.sup ){
y = plot->Y.sup;
nyh++;
}
if ( !plot->Y.force ) continue;
}
if( np==plot->maxpoints && np<n ) {
printf("reached maxpoints at %g\n", s->x[i]);
}
xp.push_back( x );
yp.push_back( y );
np = xp.size();
if( np > plot->maxpoints ) // continue in large steps
i += plot->maxpoints;
}
if ( np==0 ) {
printf( "file %zu spec %zu:"
" all %zu points outside plot range:\n",
k, j, n );
if( nxl )
printf( " %zu points < xmin\n", nxl );
if( nxh )
printf( " %zu points > xmax\n", nxh );
if( nyl )
printf( " %zu points < ymin\n", nyl );
if( nyh )
printf( " %zu points > ymax\n", nyh );
} else {
plot->plotScan( false, pstyle++,
xp, yp, e->Z(j),
e->P()->xco.str(), e->P()->yco.str(),
"data file "+strg(k)+" scan "+strg(j) );
}
/*
} else if ( fc && false ) {
// keep old vectorial version;
// might be useful for accelerated plotting
int npts = 987;
vector<double> xp(npts), yp;
for( size_t i=0; i<npts; ++i )
xp[i] = plot->X.plotcoord2value( i/(npts-1.0) );
fc->T->tree_curve_val_vec( &yp, xp, k, j );
plot->plotScan( true, cstyle++, xp, yp, e->Z(j),
e->P()->xco.str(), e->P()->yco.str(),
"curve file "+strg(k)+" scan "+strg(j) );
*/
} else if ( fc && fc->kconv!=-1 ) {
// plot only at x points of convolutand
uint kconv = fc->kconv;
uint jconv = CTree::setConvJ( k, j, kconv );
COld *fconv = NOlm::MOM[kconv].D();
if ( !fconv )
throw string( "convolutand is not a data file" );
vector<double> xp, yp;
for( size_t i=0; i<fconv->VS[jconv].size(); ++i ){
double x = fconv->VS[jconv].x[i];
if( plot->X.contained( x ) )
xp.push_back( x );
}
fc->T->tree_curve_val_vec( &yp, xp, k, j );
vector<double> xo, yo;
for( size_t i=0; i<xp.size(); ++i ){
if( plot->X.contained( xp[i] ) &&
plot->Y.contained( yp[i] ) )
continue;
xo.push_back( xp[i] );
yo.push_back( yp[i] );
}
if( xp.size()==0 ){
cout << "curve k="<<strg(k)<<", j="<<strg(j)<<" is empty\n";
return;
}
plot->plotScan( true, cstyle++, xp, yp, e->Z(j),
e->P()->xco.str(), e->P()->yco.str(), "curve file "+strg(k)+" scan "+strg(j) );
} else if ( fc ) {
// improved line plotting:
// - stop/start line when leaving/reentering the frame
// - interpolate when crossing the horizontal border (TODO)
// - interpolate at discontinuities (TODO)
double x, y;
int npts = 100;
double delta0 = 1.0/(npts-1.e-10);
double delta = delta0;
vector<CScan> SS;
CScan *s;
double pos = 0;
bool was_inside = false;
bool is_inside;
while(1){
x = plot->X.plotcoord2value( pos );
y = fc->T->tree_curve_val_sca( x, k, j );
is_inside = plot->Y.contained( y );
if( !is_inside ){
if( was_inside ){
}
goto incr_pos;
}
if( !was_inside ){
// start new line segment
SS.push_back( CScan() );
s = &( SS[SS.size()-1] );
}
s->push_back( x, y );
incr_pos:
was_inside = is_inside;
if( pos>= 1 )
break;
pos += delta;
if( pos>1 )
pos = 1;
}
for( int iS=0; iS<SS.size(); ++iS ){
plot->plotScan( true, cstyle, SS[iS].x, SS[iS].y, e->Z(j),
e->P()->xco.str(), e->P()->yco.str(),
"curve file "+strg(k)+" scan "+strg(j)+
" segment "+strg(iS) );
}
cstyle++;
} else
throw string( "PROGRAM ERROR plot invalid *e" );
}
plot->plotDoclines( e->P()->lDoc );
if( fc )
plot->plotDoclines( fc->pInfo() );
}
}