From 45a8ff9e89141ab6eb2e29cba2efadc2bedbb247 Mon Sep 17 00:00:00 2001 From: "Joachim Wuttke (laptop)" <j.wuttke@fz-juelich.de> Date: Sun, 9 May 2010 23:11:30 +0200 Subject: [PATCH] 'cr': restrict fit range several corrections to last functionality additions sorted TODO --- TODO | 125 ++++++++++++++++++------------------------- pub/CHANGELOG | 2 + pub/src/curve.cpp | 62 ++++++++++++++++++--- pub/src/curve.h | 1 + pub/src/dualplot.cpp | 17 ++++-- pub/src/dualplot.h | 1 + pub/src/file_out.cpp | 14 +++-- pub/src/frida2.cpp | 21 +++++--- pub/src/manip.cpp | 1 + pub/src/olf.h | 9 +++- 10 files changed, 155 insertions(+), 98 deletions(-) diff --git a/TODO b/TODO index 0e90102f..4a0284b1 100644 --- a/TODO +++ b/TODO @@ -1,100 +1,79 @@ == BUGS == -- as-on-disk indicator broken +- coord name algebra +- convolutand must not be defined over full energy range -== WISHLIST == +== TO CHECK == -= eilt == +- as-on-disk indicator repaired +- review definition of cauchy(,) -- coq -- Fehlerbalken in opr -- plot errorbars on/off (ge) +== WISHLIST: TODO == -WISHLIST TOFTOF -- Simultanfit -- DOS -- detailed balance -- rtof: Direkteinlese +- user interface + - restore ask callback help for lists (prompt for list, not for string) + - restore help on "?" (e.g. expression help for 'md') + - 'hf' add explanation (formulae) + - 'gp!' <filename> ? + - 'op' to give expressions for all parameters + - 'dp' output for #spec>1 is obfuscated -WISHLIST JWu -- distribute wupscat, wupsbb -- restore ask callback help for lists (prompt for list, not for string) -- restore help on "?" -- wenn oio mehrere einzelne Spektren liefert, automatisch in einen File packen -- review definition of cauchy(,) -- curve: restricted fit range -- convolution with fits ? splines ? -- convolutand must not be defined over full energy range -- mfj remove redundant doc lines -- dp output for #spec>1 is obfuscated -- coord name algebra -- pipes to avoid unnecessary files -- Debye function +- docu lines and num pars + - 'ed' should allow editing + - 'mfj' remove redundant doc lines + - 'dn' inconsistent ? + - check for completeness + - in graph file: format should allow easy read-in + - 'cp' zeigt manchmal zuviel: Anzeige ohne Par-List - customization - frida.ini location via -D - private ini file - - additional settings: graph mode; debug mode (e.g. gnuplot monitor) + - additional settings: + - let user choose points/line + - let user choose color/symbol/linestyle + - debug mode (e.g. gnuplot monitor) -- distinguish opr-functions (for modifications, with error propagation) and - and fit functions (allowing for ad-hoc additions) +- allow ad-hoc addition of fit functions + (plugin ? distinguish from opr-functions ?) + - plot accessible Q-E-space + - Debye function -SEBASTIAN +- improve existing methods + - 'oio' wenn mehrere Files mit rank=2, automatisch in einen File packen + - 'oixy': auch 'oixyd' +- distribution + - include wupscat, wupsbb -= vorrangig = -- what is dn? -- set fit range -- default empty mydefs.ini file, not overwritten by new frida version +== PROJECTS == -= wünschenswert = +- detailed balance -- legend in plot window -- oio FWHM # ? -- 2D plots using pm3d # was ist das ? -- plot accessible Q-E-space # Konzept ? -- avg, int, sum in oi -- "op" to give expressions for all parameters -- let user choose points/line -- let user choose color/symbol/linestyle -- gp <filename> # gp: ? (allg. konzept. Problem: optionale Param.) +- Simultanfit + +- DOS -= documentation and online help = +- pipes to avoid unnecessary files -- List of functions (gauss, cauchy, ...) with the corresponding expressions -- give more help -- e.g. what does "md" expect? +- avg, int, sum in oi -== ON HOLD == +- 2d plots using pm3d # was ist das ? -= gerne, wenn's jemand tut = +- 'msb' reaktivieren -- complete read i96 support +- splines (e.g. to approximate measured resolution) -= gerne, wenn ich wüsste wie = +- convolution with function + - with spline + - with TOFTOF model + + +== WAITING FOR CLUE == - Ctrl+C to abort fit - oyo; [15,j,0] --> i, y2, ... # nicht so, aber .. ?? - dirac outside conv does not really make sense # call it "resol" ?? -- normal gnuplot keys in plot window ("L" for logarithmic y scaling, zoom with mouse, ...) - -= verstehe ich noch nicht = - -- allow "for loops" in fit functions -- cut file with many z values -- allow a w-dependent expression for the resolution -- make comments in ps file complete (was somehow not complete in the apparent elastic intensity case) -- extract these commands from a ps file and "\i" it (to change limits, for example) -- command to show only fit function (and then possibly "> file.txt") -- save projects # reicht history nicht ? -- after cc, op... plot does not work -- calculate fit explicitely for export - -= überzeugt mich nicht = - -- Access to all GSL functions # unabsehbar viele, nur nach Bedarf -- allow to start bins at >0 # dann werden Anfang und Ende asymmetr. behandelt -- ". k \n" should be equal to ". \n k \n" # Konzept ? -- instead "*.2" "*::2"? # Vorteil ?? -- "exit" as synonym for "quit" -- why does "dz" show information about x values? # why not ? Gegenvorschlag ? -- let user create comments which are shown with df # wer tippt Kommentare ? +- legend in plot window +- residual plot +- Graphik über Asymptote ? Oder besser gnuplot lernen ?? \ No newline at end of file diff --git a/pub/CHANGELOG b/pub/CHANGELOG index 28bb693f..12f1fa20 100644 --- a/pub/CHANGELOG +++ b/pub/CHANGELOG @@ -4,6 +4,8 @@ Release - error propagation now working for 'm' and 'ox', 'oy', .. operations - '\wh' meta command: write history, using GNU history - 'mpaea', 'mpaer': binning that keeps error bars below given bound + - 'ge': switch error bar plotting on/off + - 'cr': set restrictions for fit range - changes in user interface: - curve evaluation expression: now c[]() instead of f[]() - 'mdy-' removed: to suppress error bars, now use oy y+-0 diff --git a/pub/src/curve.cpp b/pub/src/curve.cpp index cc0881cc..3fc7ba42 100644 --- a/pub/src/curve.cpp +++ b/pub/src/curve.cpp @@ -154,6 +154,28 @@ void NCurveFile::ChangeExpr() } +//! Set fit range restrictions. + +void NCurveFile::ChangeRange() +{ + string expr = sask( "Fit range restriction ?" ); + + NOlm::IterateC fiter; + POlc fc; + while( (fc=fiter()) ) { + fc->range_expr = expr; + fc->lDoc.push_back( "cr " + expr ); + if( expr=="" ) { + fc->range_T = PTree(); + continue; + } + PTree T; + user_xaxparse( expr.c_str(), &T ); + fc->range_T = T; + } +} + + //! Display parameter list for set of online curves. void NCurveFile::ShowPar(void) @@ -386,11 +408,11 @@ void NCurveFile::SetProperties( string which ) NOlm::IterateC fiter; POlc fc; while ( (fc=fiter()) ) { - if (which=="cwc") + if (which=="wc") fc->weighing = COlc::_LIN; - else if (which=="cwl") + else if (which=="wl") fc->weighing = COlc::_LOG; - else if (which=="cwv") + else if (which=="wv") fc->weighing = COlc::_ERR; else { printf("! unknown property"); @@ -448,7 +470,8 @@ void globalEvaluate( const double* par, int m_dat, const void *data, if ( fc->weighing==COlc::_ERR && s->dy.size()!=m_dat ) throw string( "missing error bars for weighing" ); - for ( int i=0; i<m_dat; i++) { + for ( uint i=0; i<m_dat; i++) { + if( fit[i]==INFINITY ) // TEMPORARY UNTIL exceptions can be // thrown everywhere throw string( "function evaluation return INFINITY" ); @@ -569,6 +592,27 @@ void NCurveFile::Fit( bool _allow_slow_conv ) if ( fc->ProtectedSpecs.size() ) cout << "some spectra are protected\n"; fc->lDoc.push_back( "fitted to file " + strg(fc->kd) + " " + fd->name ); + + if( fc->range_T ) { + fd = POld( fd->new_old() ); + for ( uint j=0; j<fc->nJ(); j++ ) { + PSpec ein = fd->VS(j); + PSpec eout = PSpec( new CSpec ); + eout->copy_z_base( ein ); + vector<double> range( ein->size() ); + fc->range_T->tree_vec_val( &range, 0, fc->kd, j ); + bool with_dy = ein->dy.size(); + for ( uint i=0; i<ein->size(); i++ ) { + if ( range[i] ) { + if( with_dy ) + eout->push_xyd( ein->x[i], ein->y[i], ein->dy[i] ); + else + eout->push_xy( ein->x[i], ein->y[i] ); + } + } + fd->V.push_back( eout ); + } + } data.fc = fc; data.fd = fd; @@ -580,15 +624,19 @@ void NCurveFile::Fit( bool _allow_slow_conv ) control.stepbound = NFitTune::Factor; control.maxcall = NFitTune::nCall; - for (uint j=0; j<fc->nJ(); j++) { + for ( uint j=0; j<fc->nJ(); j++ ) { if ( fc->ProtectedSpecs.contains( j ) ) continue; data.j = j; if( fiter.size()>1 ) - printf("%3d", data.k); - printf("%3d", j); + printf( "%3d", data.k ); + printf( "%3d", j ); D = fd->VS(j); C = fc->VC(j); + if ( !D->size() ){ + printf( " empty range\n" ); + continue; + } // determine number of free parameters uint npfree = 0; diff --git a/pub/src/curve.h b/pub/src/curve.h index 8fdc0462..76a769e7 100644 --- a/pub/src/curve.h +++ b/pub/src/curve.h @@ -6,6 +6,7 @@ namespace NCurveFile { void CreateFitcurve(); void CreateFreecurve(); void ChangeExpr(); + void ChangeRange(); void SetFileReference( const string& which ); void SetConvTuningPars( string which ); diff --git a/pub/src/dualplot.cpp b/pub/src/dualplot.cpp index 9be26f62..bdb209ef 100644 --- a/pub/src/dualplot.cpp +++ b/pub/src/dualplot.cpp @@ -22,7 +22,8 @@ using boost::format; //! Constructor for plot window: setup for gnuplot and postscript. CPlot::CPlot( uint _iPlot, bool _logx, bool _logy ) : - iPlot( _iPlot ), X( _logx ), Y( _logy ), maxpoints(12000) + iPlot( _iPlot ), X( _logx ), Y( _logy ), + maxpoints(12000), with_errors(true) { // Start gnuplot. Use input redirection so that gnuplot receives // commands from a gp_fifo which is created here. @@ -101,6 +102,12 @@ void CPlot::setAux( string cmd ) 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; } @@ -193,7 +200,7 @@ void CPlot::addSpec( bool as_line, int style_no, if ( as_line ) gp_fnames += str( format( " with lines lt 1 lc rgb \"#%6x\"" ) % color[ style_no % mColor ] ); - else if( dyp.size() ) + else if( with_errors && dyp.size() ) gp_fnames += " with errorbars"; FILE *gp_fd; if (!(gp_fd = fopen(gp_fnam.c_str(), "w"))) @@ -204,7 +211,7 @@ void CPlot::addSpec( bool as_line, int style_no, throw "Plot::Line: x["+strg(i)+"]="+strg(xp[i])+" out of range"; if( yp[i]<Y.inf || yp[i]>Y.sup ) throw "Plot::Line: y["+strg(i)+"]="+strg(yp[i])+" out of range"; - if( dyp.size() ) + 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]); @@ -342,8 +349,8 @@ void CPlot::gp_write( string in ) //! Format ticks and tacks for postscript file. -void CPlot::ps_ticktack(uint ntack, double *tack, int ntpt, double *ticklim, - CAxis *A) +void CPlot::ps_ticktack( + uint ntack, double *tack, int ntpt, double *ticklim, CAxis *A ) { uint i; ps_accu.push_back( "[\n" ); diff --git a/pub/src/dualplot.h b/pub/src/dualplot.h index 327f7262..0dc426e5 100644 --- a/pub/src/dualplot.h +++ b/pub/src/dualplot.h @@ -6,6 +6,7 @@ class CPlot { uint iPlot; CAxis X, Y; uint maxpoints; + bool with_errors; CPlot( uint _iPlot, bool _logx, bool _logy ); diff --git a/pub/src/file_out.cpp b/pub/src/file_out.cpp index 382352cd..ccee7e19 100644 --- a/pub/src/file_out.cpp +++ b/pub/src/file_out.cpp @@ -117,9 +117,17 @@ void NFileOut::Save_y08( FILE *file, POlo f ) PSpec s = fd->VS(j); for ( uint i=0; i<f->nZ(); i++ ) fprintf( file, " z%u: %18.10g\n", i, s->z[i] ); - fprintf( file, " xy: |2 # %u entries\n", s->size() ); - for( uint i=0; i<s->size(); i++ ) - fprintf( file, " %20.15g %20.15g\n", s->x[i], s->y[i] ); + if( s->dy.size() ) { + fprintf( file, " xyd: |2 # %u entries\n", s->size() ); + for( uint i=0; i<s->size(); i++ ) + fprintf( file, " %20.15g %20.15g %20.15g\n", + s->x[i], s->y[i], s->dy[i] ); + } else { + fprintf( file, " xy: |2 # %u entries\n", s->size() ); + for( uint i=0; i<s->size(); i++ ) + fprintf( file, " %20.15g %20.15g\n", + s->x[i], s->y[i] ); + } } else { PCurve s = fc->VC(j); for ( uint i=0; i<f->nZ(); i++ ) diff --git a/pub/src/frida2.cpp b/pub/src/frida2.cpp index a79b4b06..331f54f8 100644 --- a/pub/src/frida2.cpp +++ b/pub/src/frida2.cpp @@ -145,6 +145,7 @@ int main() " cm modify function\n" " cd set data file to fit\n" " cv set file to convolute with (cv-: no convolution)\n" + " cr set fit range restriction\n" " ci extract parameter, calculate integral, ...\n" " cp print parameters\n" "Fit to data:\n" @@ -185,6 +186,9 @@ int main() } else if (cmd == "cd" || cmd == "cv" || cmd == "cv-" ) { NCurveFile::SetFileReference( cmd.substr(1) ); + } else if (cmd == "cr" ) { + NCurveFile::ChangeRange(); + } else if (cmd == "ci") { NCurveFile::IntegralProperty(); } else if (cmd == "cp") { @@ -328,21 +332,22 @@ int main() " gnp set maximum number of points\n" ; - } else if (cmd == "gx") { + } else if ( cmd == "gx" ) { string range = sask( "Set x range" ); Plots[iPlot]->X.set_from_string( range ); - } else if (cmd == "gy") { + } else if ( cmd == "gy" ) { string range = sask( "Set y range" ); Plots[iPlot]->Y.set_from_string( range ); - } else if (cmd == "ga") { + } else if ( cmd == "ga" ) { Plots[iPlot]->X.setAuto(); Plots[iPlot]->Y.setAuto(); - } else if (cmd == "gxa") { + } else if ( cmd == "gxa" ) { Plots[iPlot]->X.setAuto(); - } else if (cmd == "gya") { + } else if ( cmd == "gya" ) { Plots[iPlot]->Y.setAuto(); - } else if (cmd.substr(0,2) == "gx" || - cmd.substr(0,2) == "gy" ) { + } else if ( cmd.substr(0,2) == "gx" || + cmd.substr(0,2) == "gy" || + cmd.substr(0,2) == "ge" ) { Plots[iPlot]->setAux( cmd ); } else if ( cmd=="gnp" ) { Plots[iPlot]->maxpoints = @@ -356,7 +361,7 @@ int main() else ps_dict = ""; Plots[iPlot]->writePostscript( ps_outdir, ps_head, ps_dict ); - } else if (cmd == "gw") { + } else if ( cmd == "gw" ) { for( uint i=0; i<nPlot; ++i ){ printf( " %c %1i %s\n", i==iPlot ? '*' : ' ', i, Plots[i]->info().c_str() ); diff --git a/pub/src/manip.cpp b/pub/src/manip.cpp index d57227a3..e4462657 100644 --- a/pub/src/manip.cpp +++ b/pub/src/manip.cpp @@ -39,6 +39,7 @@ namespace NManip { //***************************************************************************// //! Restrict forthcoming operations + void NManip::ProtectSpecs( bool on ) { NOlm::SelAssert(); diff --git a/pub/src/olf.h b/pub/src/olf.h index 258a0189..6a04e1f4 100644 --- a/pub/src/olf.h +++ b/pub/src/olf.h @@ -11,6 +11,8 @@ typedef boost::shared_ptr<class COlo> POlo; typedef boost::shared_ptr<class COld> POld; typedef boost::shared_ptr<class COlc> POlc; +typedef boost::shared_ptr<class CTree> PTree; + //! Online object: virtual base class for COld, COlc. @@ -79,7 +81,7 @@ class COlc : public COlo { vector<CCoord> PCo; // online state records: - boost::shared_ptr<class CTree> T; + PTree T; vector<CList> Fixed; enum TWgt { _LIN, _LOG, _ERR } weighing; uint kd; // internal number of data file to be fitted @@ -88,7 +90,10 @@ class COlc : public COlo { double cv_epsabs; double cv_epsrel; - COlc() : weighing( _ERR ), kconv(-1) { ; }; + string range_expr; + PTree range_T; + + COlc() : weighing( _ERR ), kconv(-1), range_expr("") { ; }; // trivially duplicated functions: PCurve VC( uint j ) const; -- GitLab