diff --git a/TODO b/TODO index 5eca0397354634c1caaf960fed1256d445fece51..10cb37c43c633c932afa14fde5f6e339af84654d 100644 --- a/TODO +++ b/TODO @@ -1,3 +1,10 @@ +==== BUGS ==== + +'make test' uses wrong executable, requires 'make install' + +get rid of ad-hoc sp_XXX functions + + ==== New functionality ==== Fourier transformed (req at TOFTOF) diff --git a/pub/ftest/mpaf.f2t b/pub/ftest/mpaf.f2t new file mode 100755 index 0000000000000000000000000000000000000000..df2d09da28e649eb0745a4d3611965cb78f9e349 --- /dev/null +++ b/pub/ftest/mpaf.f2t @@ -0,0 +1,7 @@ +#!/usr/bin/env frida +fm 10 1 h +0 oy i^2 +1 mpaf 3 +exit_unless(abs(x[2,0,2]-.7)<1e-15,"mpaf_x") +exit_unless(abs(y[2,0,2]-49-2/3)<1e-14,"mpaf_y") +exit(0) diff --git a/pub/lib/curve.cpp b/pub/lib/curve.cpp index 0b27b7cb056b2f42eff69c6625160a7af529b67e..f5a8e0cf9e7183e4fbe713baf6b68f45951db409 100644 --- a/pub/lib/curve.cpp +++ b/pub/lib/curve.cpp @@ -43,7 +43,7 @@ void NCurveFile::create_fitcurve() fc->curve_query( "Curve definition" ); fc->curve_set_defaults(); - while( ROld fd = fiter.getD() ) { + while( const COld *fd = fiter.nextD() ) { POlc fout( fc->new_olc() ); vector<double> P( fout->nP, 1.0 ); fout->name = string("fit_") + fd->name; @@ -103,7 +103,7 @@ void NCurveFile::change_expr() FileIterator fiter(NSel::selC()); - while( POlc fc = fiter.getC() ) { + while( COlc *fc = fiter.nextC() ) { int nPold = fc->nP; fc->expr = ftmp->expr; fc->T = ftmp->T; @@ -146,7 +146,7 @@ void NCurveFile::set_file_reference( const string& which ) T = user_xaxparse( expr.c_str() ); } - while ( POlc fc = fiter.getC() ) { + while ( COlc *fc = fiter.nextC() ) { int k = fiter.k(); if ( which=="v-" ) fc->kconv = -1; @@ -185,7 +185,7 @@ void NCurveFile::set_par_attr( char attr ) return; RNode Tsel = user_xaxparse( sel.c_str() ); - while( POlc fc = fiter.getC() ) { + while( COlc *fc = fiter.nextC() ) { RObjVecInt vsel = Tsel->to_index_list( fc->nP ); if (!vsel) throw S("not a valid index list"); @@ -208,7 +208,7 @@ void NCurveFile::change_range() string expr = sask( "Expression restricting the fit range" ); - while( POlc fc = fiter.getC() ) { + while( COlc *fc = fiter.nextC() ) { fc->range_expr = expr; fc->lDoc.push_back( "cr " + expr ); if( expr=="" ) { @@ -226,7 +226,7 @@ void NCurveFile::set_properties( string which ) { FileIterator fiter(NSel::selC()); - while ( POlc fc = fiter.getC() ) { + while ( COlc *fc = fiter.nextC() ) { if (which=="wc") fc->weighing = COlc::_LIN; else if (which=="wl") @@ -257,7 +257,7 @@ void NCurveFile::set_properties( string which ) void NCurveFile::pars_query() { FileIterator fiter(NSel::selC()); - while( ROlc fc = fiter.getC() ) { + while( const COlc *fc = fiter.nextC() ) { for( int ip=0; ip<fc->nP; ++ip ) NOperate::Pointwise( "p" + S(ip) ); } @@ -275,7 +275,7 @@ void NCurveFile::pars_abs() RNode Tsel = user_xaxparse( sel.c_str() ); FileIterator fiter(NSel::selC()); - while( POlc fc = fiter.getC() ) { + while( COlc *fc = fiter.nextC() ) { RObjVecInt vsel = Tsel->to_index_list( fc->nP ); if (!vsel) throw S("not a valid index list"); @@ -304,7 +304,7 @@ void NCurveFile::pars_copy() RNode Tsel = user_xaxparse( sel.c_str() ); FileIterator fiter(NSel::selC()); - while( POlc fc = fiter.getC() ) { + while( COlc *fc = fiter.nextC() ) { RObjVecInt vsel = Tsel->to_index_list( fc->nP ); if (!vsel) throw S("not a valid index list"); @@ -314,7 +314,7 @@ void NCurveFile::pars_copy() CContext ctx( fiter.k(), j ); int k2 = Tk->tree_val_idx( ctx, "K" ); int j2 = Tj->tree_val_idx( ctx, "J" ); - ROlc f2 = NOlm::mem_get_C( k2 ); + const COlc *f2 = NOlm::mem_get_C( k2 ); if( j2>=f2->nJ() ) throw "no curve " + S(j2) + " in file " + S(k2); for( int ip: vsel->v ) { @@ -331,22 +331,22 @@ void NCurveFile::pars_copy() void NCurveFile::show_pars() { FileIterator fiter(NSel::selC()); - while( ROlc fc = fiter.getC() ) + while( const COlc *fc = fiter.nextC() ) fc->print_info_table(); } //! Sets resolution file reference *sv, *kv, *jv. -void NCurveFile::get_conv( RSpec *sv, int *kv, int *jv, int k, int j ) +void NCurveFile::get_conv( const CSpec **sv, int *kv, int *jv, int k, int j ) { // set *kv and *fv for given k: - ROlc fc = NOlm::mem_get_C(k); + const COlc *fc = NOlm::mem_get_C(k); *kv = fc->kconv; if( *kv==(int)-1 ) { // no convolution *jv = -1; // should not be needed return; } - ROld fv = NOlm::mem_get_D(*kv); + const COld *fv = NOlm::mem_get_D(*kv); // set *jv and *sv for given j: if ( j>=fc->nJ() ) @@ -358,5 +358,5 @@ void NCurveFile::get_conv( RSpec *sv, int *kv, int *jv, int k, int j ) else throw "cannot convolute file " + S(k) + " with resolution " + S(*kv) + ": number of spectra does not match"; - *sv = fv->VS(*jv); + *sv = fv->VS(*jv).get(); } diff --git a/pub/lib/curve.hpp b/pub/lib/curve.hpp index 596157fcd21879c5507e71d23c893d1c8b5b5111..7f29315a9fd8bc58edfac41f25cb57899b451f62 100644 --- a/pub/lib/curve.hpp +++ b/pub/lib/curve.hpp @@ -1,7 +1,7 @@ //************************************************************************************************** -//* FRIDA: fast reliable interactive data analysis -//* (C) Joachim Wuttke 1990-, v2(C++) 2001- -//* http://apps.jcns.fz-juelich.de/frida +//* FRIDA: fast reliable interactive data analysis +//* (C) Joachim Wuttke 1990-, v2(C++) 2001- +//* http://apps.jcns.fz-juelich.de/frida //************************************************************************************************** //! \file curve.hpp @@ -26,5 +26,5 @@ namespace NCurveFile { void edit_pars(); void set_par_attr( char attr ); - void get_conv( RSpec *sv, int *kv, int *jv, int k, int j ); + void get_conv( const CSpec **sv, int *kv, int *jv, int k, int j ); } diff --git a/pub/lib/edif.cpp b/pub/lib/edif.cpp index 67512298ff871fcc55d5417710044a418a372f94..4380b905849734da9aeb36f2db98c6dbdbd9a081 100644 --- a/pub/lib/edif.cpp +++ b/pub/lib/edif.cpp @@ -38,10 +38,9 @@ using boost::format; void NEdif::show_files() { FileIterator fiter( triv::iota_list( NOlm::mem_size() ), true ); - POlo f; - while( f = fiter.get() ) { - POlc fc = P2C( f ); - POld fd = P2D( f ); + while( const COlo *f = fiter.next() ) { + const COlc *fc = dynamic_cast<const COlc*>( f ); + const COld *fd = dynamic_cast<const COld*>( f ); // prepare ztext string ztext; if( f->nZ() ){ @@ -92,11 +91,11 @@ void NEdif::show_files() void NEdif::show_spectra() { FileIterator fiter(NSel::sel()); - while( POlo f = fiter.get() ) { + while( const COlo *f = fiter.next() ) { cout << "file " << fiter.k() << ":\n"; - if ( POld fd = P2D( f ) ) { + if ( const COld *fd = dynamic_cast<const COld*>( f ) ) { printf( "%3s ", "j" ); - for ( CCoord& zco: f->ZCo ) + for ( const CCoord& zco: f->ZCo ) printf(" %12.12s", zco.str_compact().c_str()); printf( " %12s (%4s) %12s\n", "x_inf", "#pts", "x_sup" ); for ( int j=0; j<f->nJ(); j++ ) { @@ -107,7 +106,7 @@ void NEdif::show_spectra() fd->VS(j)->x[0], fd->nPts(j), fd->VS(j)->x[fd->nPts(j)-1] ); } - } else if ( POlc fc = P2C( f ) ) { + } else if ( const COlc *fc = dynamic_cast<const COlc*>( f ) ) { fc->print_info_table(); } } @@ -119,7 +118,7 @@ void NEdif::show_spectra() void NEdif::show_coord() { FileIterator fiter(NSel::sel()); - while( POlo f = fiter.get() ) { + while( const COlo *f = fiter.next() ) { string out = str( format( "%3i" ) % fiter.k() ); out += " x: " + f->xco.str_compact(); out += " y: " + f->yco.str_compact(); @@ -136,7 +135,7 @@ void NEdif::show_coord() void NEdif::show_coordZ() { FileIterator fiter(NSel::sel()); - POlo f = fiter.get(); + const COlo *f = fiter.next(); // Get current coordinates: vector<CCoord> ZCo = f->ZCo;; int nz = ZCo.size(); @@ -147,7 +146,7 @@ void NEdif::show_coordZ() printf(" there are no z coordinates\n"); // Check other input files: - while ( f = fiter.get() ) { + while ( f = fiter.next() ) { int k = fiter.k(); if ( f->ZCo.size()!=nz ) throw "file " + S(k) + " has different number of z coordinates"; @@ -163,7 +162,7 @@ void NEdif::show_coordZ() void NEdif::show_numpars() { FileIterator fiter(NSel::sel()); - while ( POlo f = fiter.get() ) { + while ( const COlo *f = fiter.next() ) { cout << "# file " << fiter.k() << ":\n"; for ( int m=0; m<f->RPar.size(); m++ ) { cout << " " << f->RPar[m].Co.str_compact() << ": " << @@ -178,7 +177,7 @@ void NEdif::show_numpars() void NEdif::show_doc() { FileIterator fiter(NSel::sel()); - while( POlo f = fiter.get() ) { + while( const COlo *f = fiter.next() ) { if( fiter.size()>1 ) cout << "# file " << fiter.k() << ":\n"; for ( int i=0; i<f->lDoc.size(); i++ ) @@ -197,7 +196,7 @@ void NEdif::show_doc() void NEdif::edit_fnam() { FileIterator fiter(NSel::sel()); - while( POlo f = fiter.get() ) { + while( COlo* f = fiter.next() ) { string fnam = wask( ("Rename " + f->name).c_str(), f->name ); if ( fnam != f->name ){ f->lDoc.push_back( "ef " + fnam + " # < " + f->name ); @@ -214,11 +213,11 @@ void NEdif::edit_coord( string which ) FileIterator fiter(NSel::sel()); // check whether all files have the same coordinate: - POlo fin = fiter.get(); + COlo *fin = fiter.next(); CCoord old_co = fin->coord( genus ); bool uniform = true; bool multifile = false; - while ( fin = fiter.get() ) { + while ( fin = fiter.next() ) { multifile = true; if ( fin->coord( genus )!=old_co ) { uniform = false; @@ -229,7 +228,7 @@ void NEdif::edit_coord( string which ) // list current values: if ( multifile && !uniform ) { fiter.reset(); - while ( fin = fiter.get() ) { + while ( fin = fiter.next() ) { cout << " file " << fiter.k() << ": " << fin->coord( genus ) <<"\n"; } } @@ -238,7 +237,7 @@ void NEdif::edit_coord( string which ) uniform ? old_co.str_compact() : "f" ); fiter.reset(); if ( multifile && in=="f" ) { - while ( fin = fiter.get() ) { + while ( fin = fiter.next() ) { old_co = fin->coord( genus ); in = sask( "Coordinate "+which+" for file "+S(fiter.k())+": name(unit)", old_co.str_compact() ); @@ -249,7 +248,7 @@ void NEdif::edit_coord( string which ) } } else { CCoord new_co(in); - while ( fin = fiter.get() ) { + while ( fin = fiter.next() ) { fin->lDoc.push_back( "ec"+string(which)+" "+new_co.str_compact()+ " # old: " + old_co.str_compact() ); fin->set_coord( genus, new_co ); @@ -263,14 +262,14 @@ void NEdif::edit_coord( string which ) void NEdif::edit_doc() { FileIterator fiter(NSel::sel()); - POlo f; + COlo *f; // obtain unused tmp file name: string fname = triv::next_tmp_file( "/tmp/frida.%i" ); // transfer to tmp file: std::ofstream Fo( fname.c_str() ); - while ( f = fiter.get() ) { + while ( f = fiter.next() ) { Fo << "### FILE " << fiter.k() << " " << f->name << " ###\n"; for ( string lin: f->lDoc ) Fo << lin << "\n"; @@ -288,7 +287,7 @@ void NEdif::edit_doc() string line; while ( triv::freadln(Fi, &line) ){ if ( line.substr(0,3)=="###" ) { - f = fiter.get(); + f = fiter.next(); if ( line.substr(3,6)!=" FILE " ) throw S("invalid change in ### line: 'FILE' not found"); string aux; @@ -319,7 +318,7 @@ void NEdif::del_numpar() { int iN = zask( "Index of numeric parameter to delete" ); FileIterator fiter(NSel::sel()); - while ( POlo f = fiter.get() ) { + while ( COlo *f = fiter.next() ) { if( iN>=f->RPar.size() ) throw "file " + S( fiter.k() ) + " has no numpar n" + S(iN); f->RPar.erase( f->RPar.begin()+iN ); diff --git a/pub/lib/expr.cpp b/pub/lib/expr.cpp index 431ca37836aee1562322b04c6573c695bc83a831..d7ebc9b0fc50c4568a21360d42e2d9a9af001ba3 100644 --- a/pub/lib/expr.cpp +++ b/pub/lib/expr.cpp @@ -38,14 +38,14 @@ CContext::CContext( int _k, int _j, int _i ) { if ( k==(int)-1 ) return; // not sure whether k=-1 makes sense - ROlo f = NOlm::mem_get(k); + const COlo *f = NOlm::mem_get(k); if ( j==(int)-1 ) { if ( f->nJ()==1 ) j = 0; } else if ( j>=f->nJ() ) throw "invalid evaluation context, j="+S(j)+">=nJ="+S(f->nJ()); - ROld fd = R2D(f); + const COld *fd = dynamic_cast<const COld*>(f); if ( !fd ) { if ( i!=(int)-1 ) throw S("invalid evaluation context, i set for curve file"); diff --git a/pub/lib/fbase.cpp b/pub/lib/fbase.cpp index efec970d2a0a5d94e5fd85cdfe2527b25bc58674..88467aca615f8a55efac4f9b91bce1371f3f8baf 100644 --- a/pub/lib/fbase.cpp +++ b/pub/lib/fbase.cpp @@ -354,7 +354,7 @@ double func_rrdm( double wred, double act0, double relwidth ) { double ret = 0; for( int i=0; i<N; ++i ) { double a = act0 + da*(i-(N-1)/2); - ret += exp( -SQR((a-act0)/actsig)/2 ) * func_cauchy( wred, exp(-a) ); + ret += exp( -SQR((a-act0)/actsig)/2 ) *func_cauchy( wred, exp(-a) ); } return ret * da / sqrt(2*PI) / actsig; } @@ -399,7 +399,7 @@ double func_ornuhl( double _w, double _a, double _b ) if ( h<DBL_MIN_EXP ) return -4; // exponent underflow f = (k+b)/(SQR(k+b)+SQR(w)); - u_next = expl( h ) * f; + u_next = expl( h ) *f; S += u; if( k==0 ) continue; diff --git a/pub/lib/file_in.cpp b/pub/lib/file_in.cpp index d760a59c825fa0f9d06b1f288b4a7f3bdb32ea39..98de9f8384c1b4d46ef04b3e986835391221c0de 100644 --- a/pub/lib/file_in.cpp +++ b/pub/lib/file_in.cpp @@ -92,18 +92,16 @@ void NFileIn::Load_yda(std::ifstream& FS, string fnam) // create output file POlo fout; - POld fd; - POlc fc; bool isdata; + COlc *fc = nullptr; if ( type == "generic tabular data" ){ isdata = true; - fd = POld( new COld ); - fout = fd; + fout = POld( new COld ); } else if ( type == "frida2 curve" ){ isdata = false; - fc = POlc( new COlc ); - fout = fc; - } else + fout = POlc( new COlc ); + fc = dynamic_cast<COlc*>(fout.get()); + } else throw "Invalid file type "+type; // read history @@ -181,7 +179,7 @@ void NFileIn::Load_yda(std::ifstream& FS, string fnam) // for curve, read formula if( !isdata ) { - string expr; + string expr; if ( version==1 ) { if( !doc["Formula"] ) throw S("DEFICIENT FILE: no formula"); diff --git a/pub/lib/file_out.cpp b/pub/lib/file_out.cpp index 8bccf2ac2f2a7bddd81f13094599339af4c184fe..c244acaa9e1e19bf8504236c7c52466c21835f41 100644 --- a/pub/lib/file_out.cpp +++ b/pub/lib/file_out.cpp @@ -25,10 +25,10 @@ #include "file_out.hpp" namespace NFileOut { - void save_yda( std::ofstream& ofs, ROlo f ); - void save_y08( FILE *file, ROlo f ); - void save_csv( FILE *file, ROlo f ); - void save_tab( FILE *file, ROlo f ); + void save_yda( std::ofstream& ofs, const COlo *f ); + void save_y08( FILE *file, const COlo *f ); + void save_csv( FILE *file, const COlo *f ); + void save_tab( FILE *file, const COlo *f ); } @@ -38,7 +38,7 @@ void NFileOut::save( string fmt, bool allow_overwrite ) { // EMBEDDED_DIALOG string outfnam; FileIterator fiter(NSel::sel()); - while( POlo f = fiter.get() ) { + while( COlo *f = fiter.next() ) { // query output file name outfnam = wask( "Save as (."+fmt+")", f->name); if ( outfnam=="" ) @@ -82,10 +82,10 @@ void NFileOut::save( string fmt, bool allow_overwrite ) //! Writes a file in format y15. -void NFileOut::save_yda( std::ofstream& ofs, ROlo f ) +void NFileOut::save_yda( std::ofstream& ofs, const COlo *f ) { - ROld fd = R2D( f ); - ROlc fc = R2C( f ); + const COld *fd = dynamic_cast<const COld*>( f ); + const COlc *fc = dynamic_cast<const COlc*>( f ); YAML::Emitter out; out << YAML::BeginMap; @@ -137,7 +137,7 @@ void NFileOut::save_yda( std::ofstream& ofs, ROlo f ) for( size_t j=0; j<f->V.size(); ++j ) { out << YAML::BeginMap; out << YAML::Key << "j" << YAML::Value << j; - RSlice s = f->V[j]; + const CSlice *s = f->V[j].get(); out << YAML::Key << "z" << YAML::Value << YAML::Flow << YAML::BeginSeq; for( auto& z: s->z ) { out << YAML::BeginMap; @@ -151,13 +151,13 @@ void NFileOut::save_yda( std::ofstream& ofs, ROlo f ) } out << YAML::EndSeq; if ( fd ) { - RSpec sd = PCAST<const CSpec>(s); + const CSpec *sd = fd->VS(j).get(); out << YAML::Key << "x" << YAML::Value << YAML::Flow << sd->x; out << YAML::Key << "y" << YAML::Value << YAML::Flow << sd->y; if( sd->has_dy() ) out << YAML::Key << "dy" << YAML::Value << YAML::Flow << sd->dy; } else { - RCurve sc = PCAST<const CCurve>(s); + const CCurve *sc = fc->VC(j).get(); out << YAML::Key << "p" << YAML::Value << YAML::Flow << sc->P; out << YAML::Key << "attr" << YAML::Value << YAML::Flow << sc->ParAttr; out << YAML::Key << "fitOutcome" << YAML::Value << sc->fitOutcome; @@ -172,10 +172,10 @@ void NFileOut::save_yda( std::ofstream& ofs, ROlo f ) //! Writes a file in format y08. -void NFileOut::save_y08( FILE *file, ROlo f ) +void NFileOut::save_y08( FILE *file, const COlo *f ) { - ROld fd = R2D( f ); - ROlc fc = R2C( f ); + const COld *fd = dynamic_cast<const COld*>( f ); + const COlc *fc = dynamic_cast<const COlc*>( f ); fprintf( file, "Meta:\n" " format: frida/y08 for yaml1\n" @@ -212,7 +212,7 @@ void NFileOut::save_y08( FILE *file, ROlo f ) for( int j=0; j<f->nJ(); j++ ){ fprintf( file, " - # table %i\n", j ); - RSlice s = f->V[j]; + CSlice *s = f->V[j].get(); for ( int i=0; i<f->nZ(); i++ ) { if( RObjNum pz = PCAST<const CObjNum>(s->z[i]) ) { fprintf( file, " z%i: %18.10g\n", i, pz->to_r() ); @@ -222,7 +222,7 @@ void NFileOut::save_y08( FILE *file, ROlo f ) fprintf( file, " z%i: %s\n", i, s->z[i]->to_s().c_str() ); } if ( fd ) { - RSpec s = fd->VS(j); + const CSpec *s = fd->VS(j).get(); if( s->dy.size() ) { fprintf( file, " xyd: |2 # %i entries\n", (int)s->size() ); for( int i=0; i<s->size(); i++ ) @@ -235,7 +235,7 @@ void NFileOut::save_y08( FILE *file, ROlo f ) s->x[i], s->y[i] ); } } else { - RCurve s = fc->VC(j); + CCurve *s = fc->VC(j).get(); for ( int i=0; i<fc->nP; i++ ) fprintf( file, " p%i: %18.10g\n", i, s->P[i] ); } @@ -245,10 +245,9 @@ void NFileOut::save_y08( FILE *file, ROlo f ) //! Writes y as tab-separated table. -void NFileOut::save_csv( FILE *file, ROlo f ) +void NFileOut::save_csv( FILE *file, const COlo *f ) { - ROld fd = R2D( f ); - ROlc fc = R2C( f ); + const COld *fd = dynamic_cast<const COld*>( f ); if( !fd ) throw S("csv save not implemented for non-data files"); for( int j=0; j<f->nJ(); j++ ){ @@ -262,10 +261,9 @@ void NFileOut::save_csv( FILE *file, ROlo f ) //! Writes data as x-y columns with z header lines. -void NFileOut::save_tab( FILE *file, ROlo f ) +void NFileOut::save_tab( FILE *file, const COlo *f ) { - ROld fd = R2D( f ); - ROlc fc = R2C( f ); + const COld *fd = dynamic_cast<const COld*>( f ); if( !fd ) throw S("tab save not implemented for non-data files"); for( int j=0; j<f->nJ(); j++ ){ diff --git a/pub/lib/fit.cpp b/pub/lib/fit.cpp index 74ca982c8ae4bc4b07e6edaee28947e5099a0d5e..97908d8043265788605ae89e36c31786fd133e5d 100644 --- a/pub/lib/fit.cpp +++ b/pub/lib/fit.cpp @@ -30,8 +30,8 @@ using boost::format; bool allow_slow_conv = true; namespace NCurveFit { - string fit_one_spec( POlc fc, ROld fd, int k, int j, const lm_control_struct& control ); - void fit_global( POlc fc, ROld fd, int k, const lm_control_struct& control ); + string fit_one_spec( COlc *fc, POld fd, int k, int j, const lm_control_struct& control ); + void fit_global( COlc *fc, POld fd, int k, const lm_control_struct& control ); } void conditionally_parallel_for_with_report( @@ -95,7 +95,7 @@ void NCurveFit::set_fit_tuning_pars( string which ) //! Returns zero when residues were computed sucessfully. static int compute_residues( - double* fvec, const int mdat, const int j, const RSpec& s, const RObj& cu, const COlc::TWgt wt ) + double* fvec, const int mdat, const int j, const CSpec *s, RObj& cu, const COlc::TWgt wt ) { try { int n = s->size(); @@ -156,8 +156,8 @@ static int compute_residues( //! Data transfer between frida and lm_minimize callback functions. typedef struct { - ROld fd; - POlc fc; + const COld *fd; + COlc *fc; int k; int j; double timeout; @@ -180,8 +180,8 @@ static void fit_evaluate( const double* par, int m_dat, const void *data, mydata = (FitDatTyp*) data; if( mydata->timeout && time(nullptr) > mydata->timeout ) throw S("reached fit time out"); - POlc fc = mydata->fc; - ROld fd = mydata->fd; + COlc *fc = mydata->fc; + const COld *fd = mydata->fd; COlc::TWgt wt = fc->weighing; PCurve c = fc->VC(mydata->j); @@ -190,7 +190,7 @@ static void fit_evaluate( const double* par, int m_dat, const void *data, if ( fc->VC(mydata->j)->ParAttr[i]!='x' ) c->P[i] = par[ip++]; - RSpec s = fd->VS(mydata->j); + const CSpec *s = fd->VS(mydata->j).get(); bool want_error = wt==COlc::_VAR || wt==COlc::_VARC; RObj cu = fc->eval_curve( s->x, mydata->k, mydata->j, want_error ); @@ -213,7 +213,7 @@ static void fit_evaluate( const double* par, int m_dat, const void *data, //! Fits one spectrum. Used when there are no global parameters. Returns an info line. -string NCurveFit::fit_one_spec( POlc fc, ROld fd, int k, int j, const lm_control_struct& control ) +string NCurveFit::fit_one_spec( COlc *fc, POld fd, int k, int j, const lm_control_struct& control ) { try { PSpec D = fd->VS(j); @@ -222,7 +222,7 @@ string NCurveFit::fit_one_spec( POlc fc, ROld fd, int k, int j, const lm_control throw S("frozen, no fit"); FitDatTyp data; data.fc = fc; - data.fd = fd; + data.fd = fd.get(); data.k = k; data.j = j; int nd = D->size(); @@ -289,8 +289,8 @@ string NCurveFit::fit_one_spec( POlc fc, ROld fd, int k, int j, const lm_control //! Data transfer between frida and lm_minimize callback functions. typedef struct { - ROld fd; - POlc fc; + const COld *fd; + COlc *fc; vector<int> J2J; vector<int> Offset; int k; @@ -314,8 +314,8 @@ static void fit_evaluate_glo( const double* par, int m_dat, const void *data, FitGloDatTyp* mydata = (FitGloDatTyp*) data; if( mydata->timeout && time(nullptr) > mydata->timeout ) throw S("reached fit time out"); - POlc fc = mydata->fc; - ROld fd = mydata->fd; + COlc *fc = mydata->fc; + const COld *fd = mydata->fd; COlc::TWgt wt = fc->weighing; int npar = 0; @@ -336,16 +336,21 @@ static void fit_evaluate_glo( const double* par, int m_dat, const void *data, bool want_error = wt==COlc::_VAR || wt==COlc::_VARC; if ( !NCurveFit::verbosity && mydata->J2J.size()>1 ) { vector<std::future<int>> workers; - for ( int jj=0; jj<mydata->J2J.size(); ++jj ) + for ( int jj=0; jj<mydata->J2J.size(); ++jj ) { workers.push_back( - std::async(std::launch::async, - [&](int _jj) -> int { - return compute_residues( - fvec, mydata->Offset[_jj], mydata->J2J[_jj], - fd->VS(mydata->J2J[_jj]), - fc->eval_curve( fd->VS(mydata->J2J[_jj])->x, - mydata->k, mydata->J2J[_jj], want_error ), - wt ); }, jj ) ); + std::async( + std::launch::async, + [&](int _jj) -> int + { + RObj cu = fc->eval_curve( fd->VS(mydata->J2J[_jj])->x, + mydata->k, mydata->J2J[_jj], want_error ); + return compute_residues( + fvec, mydata->Offset[_jj], mydata->J2J[_jj], + fd->VS(mydata->J2J[_jj]).get(), cu, wt ); + }, jj + ) + ); + } for ( auto& w: workers ) { if( w.get() ) { *userbreak = -1; @@ -354,11 +359,10 @@ static void fit_evaluate_glo( const double* par, int m_dat, const void *data, } } else { for ( int jj=0; jj<mydata->J2J.size(); ++jj ) { + RObj cu = fc->eval_curve( fd->VS(mydata->J2J[jj])->x, + mydata->k, mydata->J2J[jj], want_error ); if ( compute_residues( fvec, mydata->Offset[jj], mydata->J2J[jj], - fd->VS(mydata->J2J[jj]), - fc->eval_curve( fd->VS(mydata->J2J[jj])->x, - mydata->k, mydata->J2J[jj], want_error ), - wt ) ) { + fd->VS(mydata->J2J[jj]).get(), cu, wt ) ) { *userbreak = -1; return; } @@ -379,11 +383,11 @@ static void fit_evaluate_glo( const double* par, int m_dat, const void *data, //! Fits all spectra simultaneously. Used when there are global parameters. -void NCurveFit::fit_global( POlc fc, ROld fd, int k, const lm_control_struct& control ) +void NCurveFit::fit_global( COlc *fc, POld fd, int k, const lm_control_struct& control ) { FitGloDatTyp data; data.fc = fc; - data.fd = fd; + data.fd = fd.get(); data.k = k; // Set vector J2J to contain indices of non-frozen curves. @@ -444,7 +448,7 @@ void NCurveFit::fit_global( POlc fc, ROld fd, int k, const lm_control_struct& co fc->VC(j)->fitOutcome = (double) status.outcome; // compute deviation per spectrum: - RSpec s = fd->VS(j); + const CSpec *s = fd->VS(j).get(); RObj cu = fc->eval_curve( s->x, k, j, want_error ); int n = s->size(); double *fvec = (double*) malloc( n*sizeof(double) ); @@ -505,13 +509,13 @@ void NCurveFit::fit( bool _allow_slow_conv ) NOlm::mem_backup_store(); FileIterator fiter(NSel::selC()); - while ( POlc fc = fiter.getC() ) { + while ( COlc *fc = fiter.nextC() ) { // assert validity of parameter attributes for( int j=0; j<fc->nJ(); ++j ) fc->VC(j)->assert_attr(); - ROld fd = NOlm::mem_get_D(fc->kd); + POld fd = NOlm::sp_mem_get_D(fc->kd); // TODO get rid of shared ptr if ( fc->nJ() != fd->nJ() ) throw S("data file and curve file have different number of spectra"); if ( !fc->nJ() ) @@ -523,7 +527,7 @@ void NCurveFit::fit( bool _allow_slow_conv ) // overwrite fd with data that obye the range restriction: POld fdr( fd->new_old() ); for ( int j=0; j<fc->nJ(); j++ ) { - RSpec ein = fd->VS(j); + const CSpec *ein = fd->VS(j).get(); PSpec eout = PSpec( new CSpec ); eout->z = ein->z; vector<double> range( ein->size() ); diff --git a/pub/lib/fsel.cpp b/pub/lib/fsel.cpp index 23c1521c9f4ea07be5f794342a49e866ac468e3d..dd475b812d72cf85cb5d6f503392ee0b6c775d27 100644 --- a/pub/lib/fsel.cpp +++ b/pub/lib/fsel.cpp @@ -34,7 +34,7 @@ namespace NSel { { vector<int> ret; for( int k: FSel ) - if( P2D(NOlm::mem_get(k)) ) + if( dynamic_cast<const COld*>(NOlm::mem_get(k)) ) ret.push_back( k ); return ret; } @@ -44,13 +44,13 @@ namespace NSel { { vector<int> ret; for( int k: FSel ) - if( P2C(NOlm::mem_get(k)) ) + if( dynamic_cast<const COlc*>(NOlm::mem_get(k)) ) ret.push_back( k ); return ret; } //! First selected file. - const ROlo sel_first() + const COlo* sel_first() { if( !FSel.size() ) throw S("empty file selection"); diff --git a/pub/lib/fsel.hpp b/pub/lib/fsel.hpp index 3ef3c058a6e1ed7f66e23ef0ff5767049d4ccd41..abb0f5258d2235d353117f3f684ca72dccf83de9 100644 --- a/pub/lib/fsel.hpp +++ b/pub/lib/fsel.hpp @@ -15,7 +15,7 @@ namespace NSel { const vector<int>& sel(); const vector<int> selD(); const vector<int> selC(); - const ROlo sel_first(); + const COlo* sel_first(); const int nJ_max(); // Modify file selection: diff --git a/pub/lib/jsel.cpp b/pub/lib/jsel.cpp index c3aa18be92d972cfaa0e25cfe55d7cd1f7e143f3..184f734b35dc5fd8dd2e1eefaf227f10fe82919b 100644 --- a/pub/lib/jsel.cpp +++ b/pub/lib/jsel.cpp @@ -29,7 +29,7 @@ void JSelAsk( string quest, vector<int>& jLis ) { // EMBEDDED_DIALOG FileIterator fiter(NSel::sel()); int nJmax = 0; - while ( POlo f=fiter.get() ) + while ( const COlo *f = fiter.next() ) nJmax = std::max( nJmax, f->nJ() ); if ( nJmax==1 ){ diff --git a/pub/lib/loop.cpp b/pub/lib/loop.cpp index d11bf05a1722c113b091a5043d73e0f6bc84163f..44746a0c35dd11d5b83c33a29a6cd242dbbd2bc0 100644 --- a/pub/lib/loop.cpp +++ b/pub/lib/loop.cpp @@ -27,7 +27,7 @@ FileIterator::FileIterator( const vector<int>& _Sel, bool tolerate_empty ) reset(); } -POlo FileIterator::get() +POlo FileIterator::sp_next() { if ( iV>=Sel.size() ) return POlo(); @@ -37,12 +37,12 @@ POlo FileIterator::get() if( k>=NOlm::mem_size() ) throw "FileSelection["+S(iV)+"]="+S(k)+" exceeds online memory size "+S(NOlm::mem_size()); ++iV; - return NOlm::mem_get(k); + return NOlm::sp_mem_get(k); } -POld FileIterator::getD() +POld FileIterator::sp_nextD() { - POlo f = get(); + POlo f = sp_next(); if( !f ) return POld(); if( POld ret = P2D(f) ) @@ -50,9 +50,9 @@ POld FileIterator::getD() throw S("BUG: non-D file found by D iterator"); } -POlc FileIterator::getC() +POlc FileIterator::sp_nextC() { - POlo f = get(); + POlo f = sp_next(); if( !f ) return POlc(); if( POlc ret = P2C(f) ) diff --git a/pub/lib/loop.hpp b/pub/lib/loop.hpp index a4262b5b290592a7fb6a024f23d79f12629a9d6b..6d9aef31583907b290182ea85d8bb8b20b8ba3fb 100644 --- a/pub/lib/loop.hpp +++ b/pub/lib/loop.hpp @@ -21,7 +21,10 @@ public: int iteration() const { return iV-1; } int size() const { return Sel.size(); } // return current file, and advance iterator: - POlo get(); - POld getD(); - POlc getC(); + POlo sp_next(); + POld sp_nextD(); + POlc sp_nextC(); + COlo* next() { return sp_next().get(); } + COld* nextD() { return sp_nextD().get(); } + COlc* nextC() { return sp_nextC().get(); } }; diff --git a/pub/lib/manip.cpp b/pub/lib/manip.cpp index 6ce62113e199dc371ddfe94389d005fdcd921f9f..a825f844055a9b20528d7bf238646152f4d560dd 100644 --- a/pub/lib/manip.cpp +++ b/pub/lib/manip.cpp @@ -47,7 +47,7 @@ void NManip::freeze_slices() if ( sel=="" ) return; RNode Tsel = user_xaxparse( sel.c_str() ); - while ( POlo fio = fiter.get() ) { + while ( COlo *fio = fiter.next() ) { RObjVecInt vsel = Tsel->to_index_list( fio->nJ() ); if (!vsel) throw S("not a valid index list"); @@ -62,7 +62,7 @@ void NManip::freeze_slices() void NManip::unfreeze_all_slices() { FileIterator fiter(NSel::sel()); - while ( POlo fio = fiter.get() ) { + while ( COlo *fio = fiter.next() ) { for ( int j=0; j<fio->nJ(); ++j ) fio->V[j]->frozen = false; fio->lDoc.push_back( "m/-" ); @@ -85,11 +85,11 @@ void NManip::points_select( bool sel_del ) return; RNode Tsel = user_xaxparse( sel.c_str() ); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); fout->lDoc.push_back( "mp" + string( sel_del ? "d" : "r" ) + " " + sel ); for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); PSpec sout( new CSpec ); sout->z = sin->z; RObjVecInt vsel = Tsel->to_index_list( sin->size() ); @@ -125,12 +125,12 @@ void NManip::points_rebin() return; RNode Tsel = user_xaxparse( sel.c_str() ); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); fout->lDoc.push_back("mpa " + sel); for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); RObjVecInt vsel = Tsel->to_index_list( sin->size() ); if (!vsel) throw S("not a valid index list"); @@ -160,12 +160,12 @@ void NManip::points_rebin_by_factor() if ( ng<1 ) throw S("invalid choice"); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); fout->lDoc.push_back( "mpaf " + S(ng) ); for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin->VS(j); + const CSpec* sin = fin->VS(j).get(); int nout = sin->size() / ng; if ( nout<1 ) @@ -199,11 +199,11 @@ void NManip::points_rebin_by_err( const string& subcmd ) static double errbound; errbound = dask( "Keep error below", errbound ); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); string groupinfo; for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); PSpec sout( new CSpec() ); sout->z = sin->z; if ( !sin->has_dy() ) @@ -252,13 +252,13 @@ void NManip::points_sort() if ( T->has_dummy() ) throw S("dummy argument t not allowed when operating on data"); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { int k = fiter.k(); POld fout( fin->new_old() ); fout->lDoc.push_back("mpo "+expr); for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); PSpec sout( new CSpec ); sout->z = sin->z; @@ -286,12 +286,12 @@ void NManip::points_rebin_duplicates() { FileIterator fiter(NSel::selD()); COld fout; - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); fout->lDoc.push_back("mpq # average points when x equal"); for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); PSpec sout( new CSpec ); sout->z = sin->z; @@ -333,14 +333,14 @@ void NManip::points_symmetrize() if ( expr=="" ) return; RNode T = user_xaxparse( expr.c_str() ); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { int k = fiter.k(); POld fout( fin->new_old() ); fout->lDoc.push_back("mpsym "+expr); for (int j=0; j<fin->nJ(); j++) { int n, i, il, ih, nl, nh; - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); PSpec sout( new CSpec ); sout->z = sin->z; if ( !triv::is_ascending( sin->x ) ) @@ -390,11 +390,11 @@ void NManip::points_symmetrize() void NManip::points_remove_err() { FileIterator fiter(NSel::selD()); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); fout->lDoc.push_back("mpe-"); for (int j=0; j<fin->nJ(); j++) { - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); PSpec sout( new CSpec ); sout->z = sin->z; sout->x = sin->x; @@ -443,7 +443,7 @@ void NManip::slices_select( bool sel_del ) JSelAsk( string( sel_del ? "Delete" : "Retain" ) + " which spectra", JSel ); FileIterator fiter(NSel::sel()); - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { POlo fout( fin->new_olo() ); // JSel.evaluate( 0, fin->nJ()-1 ); fout->lDoc.push_back( "ms" + string( sel_del ? "d" : "r" ) + @@ -467,11 +467,11 @@ void NManip::slices_rebin() vector<int> JSel; static string jSel = ""; JSelAsk( "Start groups at spectra", JSel ); - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { if( fin->nJ()<=1 ) throw S("Only one slice in file"); - ROld fd = R2D( fin ); - ROlc fc = R2C( fin ); + const COld *fd = dynamic_cast<const COld*>( fin ); + const COlc *fc = dynamic_cast<const COlc*>( fin ); POlo fout( fin->new_olo() ); // JSel.evaluate( 0, fin->nJ()-1 ); if ( JSel.size()<1 || JSel[0]!=0 ) @@ -584,7 +584,7 @@ void NManip::slices_merge() static string jSel = ""; JSelAsk( "Start groups at spectra", JSel ); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); // JSel.evaluate( 0, fin->nJ()-1 ); if (JSel.size()<1 || JSel[0]!=0) @@ -624,7 +624,7 @@ void NManip::slices_spawn() if (njjIn<1) return; njj = njjIn; - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { POlo fout( fin->new_olo() ); fout->lDoc.push_back( "ms* " + S(njj) ); @@ -647,7 +647,7 @@ void NManip::slices_spawn() void NManip::exchange_x_z() { FileIterator fiter(NSel::selD()); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); int izco = fin->ZCo.size() - 1; @@ -734,7 +734,7 @@ void NManip::slices_sort() if (expr=="") return; RNode T = user_xaxparse( expr.c_str() ); - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { POlo fout( fin->new_olo() ); fout->lDoc.push_back( "mso "+expr ); @@ -757,11 +757,11 @@ void NManip::slices_sort() void NManip::slices_sort_by_z() { FileIterator fiter(NSel::sel()); - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { POlo fout( fin->new_olo() ); for (int j=0; j<fin->nJ(); j++) fout->V.push_back( fin->new_slice( j ) ); - sort( fout->V.begin(), fout->V.end(), []( const RSlice& E1, const RSlice& E2 ) { + sort( fout->V.begin(), fout->V.end(), []( const PSlice& E1, const PSlice& E2 ) { if( E1->z.size()!=E2->z.size() ) throw S("BUG: inconsistent z.size in CompareZ"); int nz = E1->z.size(); @@ -788,7 +788,7 @@ void NManip::zcoords_reorder() string com; int num1=-1, num2=-1; // nonsense initialization to suppress warning - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { nzmax = std::max( (int)nzmax, (int)fin->ZCo.size() ); } if ( nzmax<=1 ) @@ -832,7 +832,7 @@ void NManip::zcoords_reorder() } fiter.reset(); - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { POlo fout( fin->new_olo() ); fout->lDoc.push_back( com ); for (int j=0; j<fin->nJ(); j++) @@ -891,7 +891,7 @@ void NManip::zcoord_delete() return; RNode Tsel = user_xaxparse( sel.c_str() ); - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { POlo fout( fin->new_olo() ); fout->lDoc.push_back( "mz- "+sel ); @@ -932,7 +932,7 @@ void NManip::zcoord_delete() void NManip::slices_break() { FileIterator fiter(NSel::sel()); - while ( ROlo fin = fiter.get() ) { + while ( const COlo *fin = fiter.next() ) { if ( !(fin->ZCo.size()) ) throw S("no z coordinate"); CCoord zco = fin->ZCo[0]; @@ -980,28 +980,28 @@ void NManip::files_merge( const string& opts ) add_zk = true; // Determine file type, which must be the same for all input files. - ROlo fin = fiter.get(); + const COlo *fin = fiter.next(); POlo fout; - if ( ROld fd = R2D(fin) ) { + if ( dynamic_cast<const COld*>(fin) ) { fout = POld( new COld() ); - while ( fin = fiter.get() ) - if ( !R2D(fin) ) + while ( fin = fiter.next() ) + if ( !dynamic_cast<const COld*>(fin) ) throw S("cannot merge data and curve file"); - } else if ( ROlc fc = R2C(fin) ) { + } else if ( dynamic_cast<const COlc*>(fin) ) { fout = POlc( new COlc() ); - while ( fin = fiter.get() ) - if ( !R2C(fin) ) + while ( fin = fiter.next() ) + if ( !dynamic_cast<const COlc*>(fin) ) throw S("cannot merge curve and data file"); } else throw S("BUG: unexpected POlo type"); // Check and set coordinates and other file-level attributes. fiter.reset(); - fin = fiter.get(); + fin = fiter.next(); fout->xco = fin->xco; fout->yco = fin->yco; fout->ZCo = fin->ZCo; - if ( ROlc fc = R2C(fin) ) { + if ( const COlc *fc = dynamic_cast<const COlc*>(fin) ) { P2C(fout)->expr = fc->expr; P2C(fout)->evaMode = fc->evaMode; P2C(fout)->scrInpMode = fc->scrInpMode; @@ -1009,7 +1009,7 @@ void NManip::files_merge( const string& opts ) P2C(fout)->PCo = fc->PCo; P2C(fout)->T = fc->T; } - while ( fin = fiter.get() ) { + while ( fin = fiter.next() ) { if ( fin->xco != fout->xco ) throw S("different x coordinates"); if ( fin->yco != fout->yco ) @@ -1019,8 +1019,8 @@ void NManip::files_merge( const string& opts ) for ( int i=0; i<fout->ZCo.size(); ++i ) if ( fin->ZCo[i] != fout->ZCo[i] ) throw "different z"+S(i)+" coordinate"; - if ( ROlc fc = R2C(fin) ) { - if ( R2C(fin)->expr != fc->expr ) + if ( const COlc *fc = dynamic_cast<const COlc*>(fin) ) { + if ( dynamic_cast<const COlc*>(fin)->expr != fc->expr ) throw S("different curve definitions"); if ( fc->evaMode != P2C(fout)->evaMode ) throw S("different curve evaluation mode"); @@ -1033,10 +1033,10 @@ void NManip::files_merge( const string& opts ) // check RPar, and prepare move to z if they differ: fiter.reset(); - vector<CParam> RPar( fiter.get()->RPar ); + vector<CParam> RPar( fiter.next()->RPar ); int nR = RPar.size(); vector<bool> RParDiff( nR, false ); - while ( (fin = fiter.get()) ) { + while ( (fin = fiter.next()) ) { if ( fin->RPar.size() != nR ) throw S("different numer of r parameters"); for ( int i=0; i<nR; ++i ) @@ -1062,8 +1062,8 @@ void NManip::files_merge( const string& opts ) // file name: string fnam; fiter.reset(); - fnam = fiter.get()->name; - while ( fin = fiter.get() ) { + fnam = fiter.next()->name; + while ( fin = fiter.next() ) { int n2=fin->name.size(); if ( n2<fnam.size() ) fnam.erase( fnam.begin()+n2, fnam.end() ); @@ -1083,11 +1083,11 @@ void NManip::files_merge( const string& opts ) // doc lines: fout->lDoc.push_back( "mfj # " + fnam + " is merger of:" ); fiter.reset(); - ROlo f1 = fiter.get(); + const COlo *f1 = fiter.next(); fout->lDoc.push_back("- " + f1->name); for ( int i=0; i<f1->lDoc.size(); ++i ) fout->lDoc.push_back(" " + f1->lDoc[i]); - while ( (fin = fiter.get()) ) { + while ( (fin = fiter.next()) ) { fout->lDoc.push_back("- " + fin->name); string line = " "; for ( int i=0; i<fin->lDoc.size(); ++i ) { @@ -1102,7 +1102,7 @@ void NManip::files_merge( const string& opts ) // now merge the Zentries, and add new z's: fiter.reset(); - while ( fin = fiter.get() ) { + while ( fin = fiter.next() ) { for ( int j=0; j<fin->nJ(); j++ ) { PSlice eout( fin->new_slice(j) ); auto zinsert=eout->z.begin(); @@ -1129,19 +1129,19 @@ void NManip::files_merge_pointwise() { FileIterator fiter(NSel::selD()); // Copy first input file: - ROld fin = fiter.getD(); + const COld *fin = fiter.nextD(); POld fout( fin->new_old() ); fout->copy_mainvec( fin ); // Merge further input files: - while ( fin = fiter.getD() ) { + while ( fin = fiter.nextD() ) { if ( fin->nJ() != fout->nJ() ) throw S("different # spectra"); if ( fin->xco != fout->xco || fin->yco != fout->yco ) throw S("different x/y coordinates"); for ( int j=0; j<fout->nJ(); ++j ) { - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); PSpec sout = fout->VS(j); if ( sin->z.size() != sout->z.size() ) throw S("different # z entries"); @@ -1179,19 +1179,19 @@ void NManip::interpolate() RNode T = user_xaxparse( expr.c_str() ); FileIterator fiter(NSel::selD()); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { int k = fiter.k(); POld fout( fin->new_old() ); fout->copy_mainvec( fin ); CContext ctx( k ); int k2 = T->tree_val_idx( ctx, "K" ); - POld f2 = NOlm::mem_get_D( k2 ); + POld f2 = NOlm::sp_mem_get_D( k2 ); fout->lDoc.push_back( "mi " + S(k2) ); for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin ->VS(j); + const CSpec *sin = fin ->VS(j).get(); PSpec s2 = f2 ->VS(j); PSpec sout= fout->VS(j); int i2 = 0; @@ -1235,9 +1235,9 @@ void NManip::make_histogram() double xmin=0, xmax=0; // initialized to prevent warning int npts = 0; FileIterator fiter(NSel::selD()); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin ->VS(j); + const CSpec *sin = fin ->VS(j).get(); for ( int i=0; i<sin->size(); i++ ) { double x = sin->x[i]; if ( !npts ) { @@ -1271,12 +1271,12 @@ void NManip::make_histogram() // Create and fill output histograms: fiter.reset(); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); fout->lDoc.push_back("mhm " + S(nbin) ); for ( int j=0; j<fin->nJ(); j++ ) { - RSpec sin = fin->VS(j); + const CSpec *sin = fin->VS(j).get(); PSpec sout( new CSpec() ); sout->z = sin->z; sout->x = xout; diff --git a/pub/lib/mem.cpp b/pub/lib/mem.cpp index 6c8a0b6a36d0cd70a2e74f45aeb8c6b35a193875..5a52173671605f9bbe243e6695cdfe3eb24ae026 100644 --- a/pub/lib/mem.cpp +++ b/pub/lib/mem.cpp @@ -42,7 +42,7 @@ namespace NOlm { } //! Return pointer to file with index k. - const POlo mem_get( int k ) + POlo sp_mem_get( int k ) { if ( k<0 ) throw "negative number " + S(k) + " cannot be a file index"; @@ -51,24 +51,30 @@ namespace NOlm { return MOM[k]; } + const COlo* mem_get( int k ) { return sp_mem_get(k).get(); } + //! Return pointer to file with index k, provided it contains data. - const POld mem_get_D( int k ) + POld sp_mem_get_D( int k ) { - POld ret = P2D( mem_get( k ) ); + POld ret = P2D( sp_mem_get( k ) ); if( !ret ) throw "file " + S(k) + " is not a data file"; return ret; } + const COld* mem_get_D( int k ) { return sp_mem_get_D(k).get(); } + //! Return pointer to file with index k, provided it contains curves. - const POlc mem_get_C( int k ) + POlc sp_mem_get_C( int k ) { - POlc ret = P2C( mem_get( k ) ); + POlc ret = P2C( sp_mem_get( k ) ); if( !ret ) throw "file " + S(k) + " is not a curve file"; return ret; } + const COlc* mem_get_C( int k ) { return sp_mem_get_C(k).get(); } + //************************************************************************************************** //* Modify online memory, and adapt file NSel::selection @@ -105,7 +111,7 @@ namespace NOlm { void mem_copy() { FileIterator fiter(NSel::sel()); - while ( ROlo fin=fiter.get() ){ + while ( const COlo *fin=fiter.next() ){ POlo fout( fin->new_olo( false ) ); fout->copy_mainvec( fin ); mem_store( fout ); @@ -129,7 +135,7 @@ namespace NOlm { { BAK.clear(); FileIterator fiter(NSel::sel()); - while ( ROlo fin=fiter.get() ){ + while ( const COlo *fin=fiter.next() ){ POlo fout( fin->new_olo( false ) ); fout->copy_mainvec( fin ); BAK.push_back( fout ); diff --git a/pub/lib/mem.hpp b/pub/lib/mem.hpp index a0bec5f378a3f4e4318672351079527e613c1383..e5ac8ac6ed9f2753048b4aefd0827ebc3cb19d22 100644 --- a/pub/lib/mem.hpp +++ b/pub/lib/mem.hpp @@ -12,9 +12,11 @@ namespace NOlm { // Read from online memory: - const POlo mem_get ( int k ); - const POld mem_get_D( int k ); - const POlc mem_get_C( int k ); + POlo sp_mem_get ( int k ); // TODO wegdamit + const COlo* mem_get ( int k ); + const COld* mem_get_D( int k ); + POld sp_mem_get_D( int k ); // TODO get rid of this workaround + const COlc* mem_get_C( int k ); const int mem_size(); // Modify online memory (file selection is changed accordingly): diff --git a/pub/lib/node.cpp b/pub/lib/node.cpp index 6c7a10c03e0776abe1129ab54963fd8867a1e031..2f1505eedbfda9d3b1552d1ae2238a0f8ee27555 100644 --- a/pub/lib/node.cpp +++ b/pub/lib/node.cpp @@ -637,8 +637,8 @@ RObj CNodeCvin::tree_val( const CContext& ctx ) const } } int k = ref->get_k( ctx ); - ROlo f = NOlm::mem_get(k); - ROlc fc = R2C( f ); + const COlo *f = NOlm::mem_get(k); + const COlc *fc = dynamic_cast<const COlc*>( f ); if ( !fc ) throw S("This is not a curve file => cannot integrate"); int j = PCAST<const CRef2>(ref)->get_j( ctx, f->nJ() ); @@ -779,25 +779,25 @@ CCoord CNodeFile::node_coord( int _k ) const } else { k = _k; } - ROlo f = NOlm::mem_get(k); + const COlo *f = NOlm::mem_get(k); return file_coord( f ); } RObj CNodeFileNJ::tree_val_scalar( const CContext& ctx ) const { - ROlo f = NOlm::mem_get( ref->get_k( ctx ) ); + const COlo *f = NOlm::mem_get( ref->get_k( ctx ) ); return RObjInt( new CObjInt( f->nJ() ) ); } RObj CNodeFileR::tree_val_scalar( const CContext& ctx ) const { - ROlo f = NOlm::mem_get( ref->get_k( ctx ) ); + const COlo *f = NOlm::mem_get( ref->get_k( ctx ) ); if( num<0 || num>=f->RPar.size() ) throw "file has no parameter r"+S(num); return RObjDbl( new CObjDbl( f->RPar[num].val ) ); } -CCoord CNodeFileR::file_coord( ROlo f ) const +CCoord CNodeFileR::file_coord( const COlo *f ) const { if( num<0 || num>=f->RPar.size() ) throw "file has no parameter r"+S(num); @@ -812,15 +812,15 @@ CCoord CNodeFileR::file_coord( ROlo f ) const RObj CNodeSliceZ::tree_val_scalar( const CContext& ctx ) const { - ROlo f = NOlm::mem_get( ref->get_k( ctx ) ); + const COlo *f = NOlm::mem_get( ref->get_k( ctx ) ); int j = PCAST<const CRef2>(ref)->get_j( ctx, f->nJ() ); - RSlice s = f->V[j]; + const CSlice * s = f->V[j].get(); if ( num<0 || num>=s->z.size() ) throw "z"+S(num)+" requested, while only "+S(s->z.size())+" z entries are defined"; return s->z[num]; } -CCoord CNodeSliceZ::file_coord( ROlo f ) const +CCoord CNodeSliceZ::file_coord( const COlo *f ) const { if ( num<0 || num>=f->ZCo.size() ) throw "z"+S(num)+" requested, while only "+S(f->ZCo.size())+" z coordinates are defined"; @@ -829,8 +829,8 @@ CCoord CNodeSliceZ::file_coord( ROlo f ) const RObj CNodeSpecNI::tree_val_scalar( const CContext& ctx ) const { - ROlo f = NOlm::mem_get( ref->get_k( ctx ) ); - ROld fd = R2D( f ); + const COlo *f = NOlm::mem_get( ref->get_k( ctx ) ); + const COld *fd = dynamic_cast<const COld*>( f ); if ( !fd ) throw f->type()+" file has no property "+tree_info(); int j = PCAST<const CRef2>(ref)->get_j( ctx, f->nJ() ); @@ -839,25 +839,25 @@ RObj CNodeSpecNI::tree_val_scalar( const CContext& ctx ) const RObj CNodeCurve::tree_val_scalar( const CContext& ctx ) const { - ROlo f = NOlm::mem_get( ref->get_k( ctx ) ); - ROlc fc = R2C( f ); + const COlo *f = NOlm::mem_get( ref->get_k( ctx ) ); + const COlc *fc = dynamic_cast<const COlc*>( f ); if ( !fc ) throw f->type()+" file has no property "+tree_info(); int j = PCAST<const CRef2>(ref)->get_j( ctx, f->nJ() ); - return curve_val_scalar( fc->VC(j) ); + return curve_val_scalar( fc->VC(j).get() ); } -RObj CNodeCurveP::curve_val_scalar( RCurve c ) const +RObj CNodeCurveP::curve_val_scalar( const CCurve * c ) const { if ( num<0 || num>=c->P.size() ) throw "curve has no p"+S(num); return RObjDbl( new CObjDbl( c->P[num] ) ); } -CCoord CNodeCurveP::file_coord( ROlo f ) const +CCoord CNodeCurveP::file_coord( const COlo *f ) const { - ROlc fc = R2C(f); + const COlc *fc = dynamic_cast<const COlc*>(f); if ( !fc ) throw f->type()+" file has no property "+tree_info(); if ( num<0 || num>=fc->PCo.size() ) @@ -866,17 +866,17 @@ CCoord CNodeCurveP::file_coord( ROlo f ) const } -RObj CNodeCurveOutcome::curve_val_scalar( RCurve c ) const +RObj CNodeCurveOutcome::curve_val_scalar( const CCurve * c ) const { return RObjInt( new CObjInt( c->fitOutcome ) ); } -RObj CNodeCurveChi2::curve_val_scalar( RCurve c ) const +RObj CNodeCurveChi2::curve_val_scalar( const CCurve * c ) const { return RObjDbl( new CObjDbl( c->fitChi2 ) ); } -RObj CNodeCurveR2::curve_val_scalar( RCurve c ) const +RObj CNodeCurveR2::curve_val_scalar( const CCurve * c ) const { return RObjDbl( new CObjDbl( c->fitR2 ) ); } @@ -916,12 +916,12 @@ RObj CNodePoint::tree_val( const CContext& ctx ) const RObj CNodePoint::tree_val_point( const CContext& ctx, int i ) const { - ROlo f = NOlm::mem_get( ref->get_k( ctx ) ); - ROld fd = R2D( f ); + const COlo *f = NOlm::mem_get( ref->get_k( ctx ) ); + const COld *fd = dynamic_cast<const COld*>( f ); if ( !fd ) throw f->type()+" file has no property "+tree_info(); int j = PCAST<const CRef2>(ref)->get_j( ctx, f->nJ() ); - RSpec s = fd->VS(j); + const CSpec *s = fd->VS(j).get(); if ( i<0 || i>= s->size() ) throw "invalid data point index "+S(i)+", spectrum has length "+S(s->size()); return spec_val_point( s, i, ctx.want_error ); @@ -932,12 +932,12 @@ RObj CNodePoint::tree_val_vector( const CContext& ctx ) const { if ( ctx.dim!=CContext::_VI ) throw S("BUG: unexpected context in tree_val_vector"); - ROlo f = NOlm::mem_get( ref->get_k( ctx ) ); - ROld fd = R2D( f ); + const COlo *f = NOlm::mem_get( ref->get_k( ctx ) ); + const COld *fd = dynamic_cast<const COld*>( f ); if ( !fd ) throw f->type()+" file has no property "+tree_info(); int j = PCAST<const CRef2>(ref)->get_j( ctx, f->nJ() ); - RSpec s = fd->VS(j); + const CSpec *s = fd->VS(j).get(); int n = ctx.nv; if ( n<0 ) // ctx.nv=-1 is used in integral operations n = s->size(); @@ -947,33 +947,33 @@ RObj CNodePoint::tree_val_vector( const CContext& ctx ) const } -RObj CNodePointX::spec_val_point( RSpec s, int i, bool want_error ) const +RObj CNodePointX::spec_val_point( const CSpec *s, int i, bool want_error ) const { return RObjDbl( new CObjDbl( s->x[i] ) ); } -RObj CNodePointX::spec_val_vector( RSpec s, int n, bool want_error ) const +RObj CNodePointX::spec_val_vector( const CSpec *s, int n, bool want_error ) const { return RObjVecDbl( new CObjVecDbl( vector<double>(s->x.begin(), s->x.begin()+n ) ) ); } -CCoord CNodePointX::file_coord( ROlo f ) const +CCoord CNodePointX::file_coord( const COlo *f ) const { - ROld fd = R2D( f ); + const COld *fd = dynamic_cast<const COld*>( f ); if ( !fd ) throw f->type()+" file has no property "+tree_info(); return fd->xco; } -RObj CNodePointY::spec_val_point( RSpec s, int i, bool want_error ) const +RObj CNodePointY::spec_val_point( const CSpec *s, int i, bool want_error ) const { return want_error && s->has_dy() ? RObjEnu( new CObjEnu( s->y[i], s->dy[i] ) ) : RObjDbl( new CObjDbl( s->y[i] ) ); } -RObj CNodePointY::spec_val_vector( RSpec s, int n, bool want_error ) const +RObj CNodePointY::spec_val_vector( const CSpec *s, int n, bool want_error ) const { return want_error && s->has_dy() ? RObjVecEnu( new CObjVecEnu( vector<double>(s->y.begin(), s->y.begin()+n ), @@ -981,30 +981,30 @@ RObj CNodePointY::spec_val_vector( RSpec s, int n, bool want_error ) const RObjVecDbl( new CObjVecDbl( vector<double>(s->y.begin(), s->y.begin()+n ) ) ); } -CCoord CNodePointY::file_coord( ROlo f ) const +CCoord CNodePointY::file_coord( const COlo *f ) const { - ROld fd = R2D( f ); + const COld *fd = dynamic_cast<const COld*>( f ); if ( !fd ) throw f->type()+" file has no property "+tree_info(); return fd->yco; } -RObj CNodePointDY::spec_val_point( RSpec s, int i, bool want_error ) const +RObj CNodePointDY::spec_val_point( const CSpec *s, int i, bool want_error ) const { return RObjDbl( new CObjDbl( s->dy[i] ) ); } -RObj CNodePointDY::spec_val_vector( RSpec s, int n, bool want_error ) const +RObj CNodePointDY::spec_val_vector( const CSpec *s, int n, bool want_error ) const { if ( !s->has_dy() ) throw S("requested dy does not exist"); return RObjVecDbl( new CObjVecDbl( vector<double>(s->dy.begin(), s->dy.begin()+n ) ) ); } -CCoord CNodePointDY::file_coord( ROlo f ) const +CCoord CNodePointDY::file_coord( const COlo *f ) const { - ROld fd = R2D( f ); + const COld *fd = dynamic_cast<const COld*>( f ); if ( !fd ) throw f->type()+" file has no property "+tree_info(); return CCoord( "d"+fd->yco.name, fd->yco.unit ); @@ -1056,7 +1056,7 @@ RObj CNodeCev::tree_val( const CContext& ctx ) const try { // get k, curve, j from context: int k = ref->get_k( ctx ); - ROlc fc = NOlm::mem_get_C( k ); + const COlc *fc = NOlm::mem_get_C( k ); int j = PCAST<const CRef2>(ref)->get_j( ctx, fc->nJ() ); // evaluate function argument: @@ -1093,7 +1093,7 @@ RObj CNodeMixin::tree_val( const CContext& ctx ) const // Resolution is given by context 'ctx': int kconv, jconv; - RSpec sv; + const CSpec *sv; NCurveFile::get_conv( &sv, &kconv, &jconv, ctx.k, ctx.j ); PObjVecDbl pret( ctx.want_error ? (new CObjVecDbl) : (new CObjVecEnu) ); @@ -1122,7 +1122,7 @@ RObj CNodeConv::copy_theory( const CContext& ctx, double theshift ) const return theory->tree_val( ctx ); } -RObjVecDbl CNodeConv::convolve( const CContext& ctx, double theshift, const RSpec& sv, +RObjVecDbl CNodeConv::convolve( const CContext& ctx, double theshift, const CSpec* sv, double conv_norm, double conv_step ) const { int n = ctx.vt->size(); @@ -1206,7 +1206,7 @@ RObj CNodePConv::copy_theory( const CContext& ctx, double theshift ) const // for which we should go to shared pointers ... } -RObjVecDbl CNodePConv::convolve( const CContext& ctx, double theshift, const RSpec& sv, +RObjVecDbl CNodePConv::convolve( const CContext& ctx, double theshift, const CSpec* sv, double conv_norm, double conv_step ) const { if( conv_norm <= 0 ) @@ -1333,7 +1333,7 @@ RObj CNodeDirac::copy_theory( const CContext& ctx, double theshift ) const //! Convolution or Dirac function: conv(f, shift), dirac(shift). -RObjVecDbl CNodeDirac::convolve( const CContext& ctx, double theshift, const RSpec& sv, +RObjVecDbl CNodeDirac::convolve( const CContext& ctx, double theshift, const CSpec* sv, double conv_norm, double conv_step ) const { int n = ctx.vt->size(); diff --git a/pub/lib/node.hpp b/pub/lib/node.hpp index e338c79092b3871d77242f2c99ca12fabe4464a2..4a49a812780c166365ffe4c0a36a1d578498f82f 100644 --- a/pub/lib/node.hpp +++ b/pub/lib/node.hpp @@ -203,7 +203,7 @@ class CNodeFile: public CNode { private: virtual RObj tree_val_scalar( const CContext& ctx ) const { throw S("BUG: unforeseen call to CNodeFile::tree_val_scalar"); } - virtual CCoord file_coord( ROlo f ) const { return CCoord(name(), ""); } + virtual CCoord file_coord( const COlo *f ) const { return CCoord(name(), ""); } protected: RRef ref; public: @@ -234,7 +234,7 @@ class CNodeFileR: public CNodeFile { int num; public: CNodeFileR( int _num, RRef _ref ) : CNodeFile( _ref), num(_num) {} - CCoord file_coord( ROlo f ) const; + CCoord file_coord( const COlo *f ) const; string name() const { return "r" + S(num); } RObj tree_val_scalar( const CContext& ctx ) const; }; @@ -255,7 +255,7 @@ class CNodeSliceZ: public CNodeSlice { int num; public: CNodeSliceZ( int _num, RRef _ref ) : CNodeSlice( _ref), num(_num) {} - CCoord file_coord( ROlo f ) const; + CCoord file_coord( const COlo *f ) const; string name() const { return "z" + S(num); } RObj tree_val_scalar( const CContext& ctx ) const; }; @@ -277,7 +277,7 @@ class CNodeCurve: public CNodeSlice { public: CNodeCurve( RRef _ref ) : CNodeSlice( _ref ) {} RObj tree_val_scalar( const CContext& ctx ) const; - virtual RObj curve_val_scalar( RCurve c ) const = 0; + virtual RObj curve_val_scalar( const CCurve * c ) const = 0; }; @@ -288,9 +288,9 @@ class CNodeCurveP: public CNodeCurve { int num; public: CNodeCurveP( int _num, RRef _ref ) : CNodeCurve( _ref), num(_num) {} - CCoord file_coord( ROlo f ) const; + CCoord file_coord( const COlo *f ) const; string name() const { return "p" + S(num); } - RObj curve_val_scalar( RCurve c ) const; + RObj curve_val_scalar( const CCurve * c ) const; int npar() const { return num+1; } }; @@ -301,7 +301,7 @@ class CNodeCurveOutcome: public CNodeCurve { public: CNodeCurveOutcome( RRef _ref ) : CNodeCurve( _ref) {} string name() const { return "oc"; } - RObj curve_val_scalar( RCurve c ) const; + RObj curve_val_scalar( const CCurve * c ) const; }; @@ -311,7 +311,7 @@ class CNodeCurveChi2: public CNodeCurve { public: CNodeCurveChi2( RRef _ref ) : CNodeCurve( _ref) {} string name() const { return "chi2"; } - RObj curve_val_scalar( RCurve c ) const; + RObj curve_val_scalar( const CCurve * c ) const; }; @@ -321,7 +321,7 @@ class CNodeCurveR2: public CNodeCurve { public: CNodeCurveR2( RRef _ref ) : CNodeCurve( _ref) {} string name() const { return "R2"; } - RObj curve_val_scalar( RCurve c ) const; + RObj curve_val_scalar( const CCurve * c ) const; }; @@ -332,8 +332,8 @@ class CNodePoint: public CNodeSlice { RObj tree_val_point( const CContext& ctx, int i ) const; RObj tree_val_vector( const CContext& ctx ) const; protected: - virtual RObj spec_val_point( RSpec s, int i, bool want_error ) const = 0; - virtual RObj spec_val_vector( RSpec s, int n, bool want_error ) const = 0; + virtual RObj spec_val_point( const CSpec *s, int i, bool want_error ) const = 0; + virtual RObj spec_val_vector( const CSpec *s, int n, bool want_error ) const = 0; public: CNodePoint( RRef _ref ) : CNodeSlice( _ref ) {} RObj tree_val( const CContext& ctx ) const; @@ -345,10 +345,10 @@ class CNodePoint: public CNodeSlice { class CNodePointX: public CNodePoint { public: CNodePointX( RRef _ref ) : CNodePoint( _ref ) {} - CCoord file_coord( ROlo f ) const; + CCoord file_coord( const COlo *f ) const; string name() const { return "x"; } - RObj spec_val_point( RSpec s, int i, bool want_error ) const; - RObj spec_val_vector( RSpec s, int n, bool want_error ) const; + RObj spec_val_point( const CSpec *s, int i, bool want_error ) const; + RObj spec_val_vector( const CSpec *s, int n, bool want_error ) const; }; @@ -357,10 +357,10 @@ class CNodePointX: public CNodePoint { class CNodePointY: public CNodePoint { public: CNodePointY( RRef _ref ) : CNodePoint( _ref ) {} - CCoord file_coord( ROlo f ) const; + CCoord file_coord( const COlo *f ) const; string name() const { return "y"; } - RObj spec_val_point( RSpec s, int i, bool want_error ) const; - RObj spec_val_vector( RSpec s, int n, bool want_error ) const; + RObj spec_val_point( const CSpec *s, int i, bool want_error ) const; + RObj spec_val_vector( const CSpec *s, int n, bool want_error ) const; }; @@ -369,10 +369,10 @@ class CNodePointY: public CNodePoint { class CNodePointDY: public CNodePoint { public: CNodePointDY( RRef _ref ) : CNodePoint( _ref ) {} - CCoord file_coord( ROlo f ) const; + CCoord file_coord( const COlo *f ) const; string name() const { return "dy"; } - RObj spec_val_point( RSpec s, int i, bool want_error ) const; - RObj spec_val_vector( RSpec s, int n, bool want_error ) const; + RObj spec_val_point( const CSpec *s, int i, bool want_error ) const; + RObj spec_val_vector( const CSpec *s, int n, bool want_error ) const; }; @@ -432,7 +432,7 @@ class CNodeMixin: public CNode { bool i_dependent() const { return shift->i_dependent(); } RObj tree_val( const CContext& ctx ) const; virtual RObj copy_theory( const CContext& ctx, double theshift ) const =0; - virtual RObjVecDbl convolve( const CContext& ctx, double theshift, const RSpec& sv, + virtual RObjVecDbl convolve( const CContext& ctx, double theshift, const CSpec* sv, double conv_norm, double conv_step ) const =0; bool has_conv() const { return true; } }; @@ -461,7 +461,7 @@ class CNodeConv: public CNodeConvBase { CNodeConv( const RNode& _theory, const RNode& _shift=RNode( new CNodeVal( 0.0 ) ) ) : CNodeConvBase( _theory, _shift ) {} RObj copy_theory( const CContext& ctx, double theshift ) const; - RObjVecDbl convolve( const CContext& ctx, double theshift, const RSpec& sv, + RObjVecDbl convolve( const CContext& ctx, double theshift, const CSpec *sv, double conv_norm, double conv_step ) const; string tree_info() const { return "conv(" + theory->tree_info() +","+ shift->tree_info() + ")"; } @@ -475,7 +475,7 @@ class CNodePConv: public CNodeConvBase { CNodePConv( const RNode& _theory, const RNode& _shift=RNode( new CNodeVal( 0.0 ) ) ) : CNodeConvBase( _theory, _shift ) {}; RObj copy_theory( const CContext& ctx, double theshift ) const; - RObjVecDbl convolve( const CContext& ctx, double theshift, const RSpec& sv, + RObjVecDbl convolve( const CContext& ctx, double theshift, const CSpec* sv, double conv_norm, double conv_step ) const; string tree_info() const { return "pconv(" + theory->tree_info() +","+ shift->tree_info() + ")"; } @@ -489,7 +489,7 @@ class CNodeDirac: public CNodeMixin { CNodeDirac( const RNode& _shift=RNode( new CNodeVal( 0.0 ) ) ) : CNodeMixin( _shift ) {}; RObj copy_theory( const CContext& ctx, double theshift ) const; - RObjVecDbl convolve( const CContext& ctx, double theshift, const RSpec& sv, + RObjVecDbl convolve( const CContext& ctx, double theshift, const CSpec *sv, double conv_norm, double conv_step ) const; bool has_dummy() const { return true; } // has implicit t dependence CCoord node_coord( int k ) const; diff --git a/pub/lib/olf.cpp b/pub/lib/olf.cpp index 5b5e037d94a0c646b1bb06aab01950de4bb1ae20..146816fd13e9da756e0d0351db1336ed258d8de4 100644 --- a/pub/lib/olf.cpp +++ b/pub/lib/olf.cpp @@ -104,7 +104,7 @@ PSlice COlc::new_slice( int j ) const POlo COld::new_olo( bool modified ) const { POlo ret( new COld( *this ) ); - ret->V.clear(); + ret->V.clear(); // TODO: this is inefficient; don't copy then delete. ret->as_on_disk = as_on_disk && !modified; return ret; } @@ -130,9 +130,9 @@ POlc COlc::new_olc( bool modified ) const //! Set main vector by copying from fin. -void COld::copy_mainvec( ROlo fin ) +void COld::copy_mainvec( const COlo *fin ) { - ROld fd = R2D( fin ); + const COld *fd = dynamic_cast<const COld*>( fin ); V.resize( fin->nJ() ); for( int j=0; j<fin->nJ(); ++j ) V[j] = PSlice( new CSpec( *(fd->VS(j)) ) ); @@ -141,9 +141,9 @@ void COld::copy_mainvec( ROlo fin ) //! Set main vector by copying from fin. -void COlc::copy_mainvec( ROlo fin ) +void COlc::copy_mainvec( const COlo *fin ) { - ROlc fc = R2C( fin ); + const COlc *fc = dynamic_cast<const COlc*>( fin ); V.resize( fin->nJ() ); for( int j=0; j<fin->nJ(); ++j ) V[j] = PSlice( new CCurve( *(fc->VC(j)) ) ); @@ -332,7 +332,7 @@ int COld::nPts() const bool COld::has_nonzero_dy() const { return std::all_of( V.cbegin(), V.cend(), - [](RSlice s){ return (PCAST<const CSpec>(s))->has_nonzero_dy(); } ); + [](const PSlice& s){ return (PCAST<const CSpec>(s))->has_nonzero_dy(); } ); } //************************************************************************************************** diff --git a/pub/lib/olf.hpp b/pub/lib/olf.hpp index b03ba567ac7beafa18cee26da198a7cf37c77498..f57777cf3821706d7db2e4679f47f21761fa9777 100644 --- a/pub/lib/olf.hpp +++ b/pub/lib/olf.hpp @@ -34,7 +34,7 @@ class COlo { string info_line( int j ) const; void check_integrity() const; - virtual void copy_mainvec( ROlo fin ) =0; + virtual void copy_mainvec( const COlo *fin ) =0; int nJ() const { return V.size(); } int nZ() const { return ZCo.size(); } @@ -59,7 +59,7 @@ class COld : public COlo { PSlice new_slice( int j ) const; POlo new_olo( bool modified=true ) const; POld new_old( bool modified=true ) const; - void copy_mainvec( ROlo fin ); + void copy_mainvec( const COlo *fin ); // Overloaded functions: void set_coord( const CGenus& genus, CCoord& co ); @@ -114,7 +114,7 @@ public: PSlice new_slice( int j ) const; POlo new_olo( bool modified=true ) const; POlc new_olc( bool modified=true ) const; - void copy_mainvec( ROlo fin ); + void copy_mainvec( const COlo *fin ); // Overloaded functions: void set_coord( const CGenus& genus, CCoord& co ); @@ -141,6 +141,4 @@ public: //! Casts. #define P2D(p) std::dynamic_pointer_cast< class COld>(p) -#define R2D(p) std::dynamic_pointer_cast<const class COld>(p) #define P2C(p) std::dynamic_pointer_cast< class COlc>(p) -#define R2C(p) std::dynamic_pointer_cast<const class COlc>(p) diff --git a/pub/lib/opr.cpp b/pub/lib/opr.cpp index bc378d583c5fec5c257d301804ef22747adf75fc..fb993ae598990fd72c5d20a2228161fe2bd1c908 100644 --- a/pub/lib/opr.cpp +++ b/pub/lib/opr.cpp @@ -48,7 +48,7 @@ void NOperate::show( const string& subcmd ) if( T->has_dummy() ) throw S("dummy argument t not allowed when operating on data"); - while( ROld fin = fiter.getD() ) { + while( const COld *fin = fiter.nextD() ) { if ( fiter.size()>1 ) printf( ".. file %i:\n", fiter.k() ); @@ -65,7 +65,7 @@ void NOperate::show( const string& subcmd ) outp = ""; for( int j=0; j<fin->nJ(); j++ ){ - RSpec s = fin->VS(j); + const CSpec *s = fin->VS(j).get(); if( i>= s->size() ) continue; if ( ! T->tree_point_dbl( fiter.k(), j, i ) ) continue; @@ -106,16 +106,16 @@ void NOperate::select( bool sel_del ) if( T->has_dummy() ) throw S("dummy argument t not allowed when operating on data"); - while ( ROlo fin=fiter.get() ) { - ROld fd = R2D( fin ); - ROlc fc = R2C( fin ); + while ( const COlo *fin=fiter.next() ) { + const COld *fd = dynamic_cast<const COld*>( fin ); + const COlc *fc = dynamic_cast<const COlc*>( fin ); POlo fout( fin->new_olo() ); fout->lDoc.push_back((sel_del ? "md " : "mr ") + expr); for ( int j=0; j<fin->nJ(); j++ ) { if ( fd ) { - RSpec ein = fd->VS(j); + const CSpec *ein = fd->VS(j).get(); PSpec eout = PSpec( new CSpec ); eout->z = ein->z; vector<double> lvec( ein->size() ); @@ -160,10 +160,10 @@ void NOperate::Pointwise( string llabel ) if( T->has_dummy() ) throw S("dummy argument t not allowed when operating on data"); - while( ROlo fin=fiter.get() ){ + while( const COlo *fin=fiter.next() ){ int k = fiter.k(); - ROld fd = R2D( fin ); - ROlc fc = R2C( fin ); + const COld *fd = dynamic_cast<const COld*>( fin ); + const COlc *fc = dynamic_cast<const COlc*>( fin ); if ( lref.pointwise() && !fd ) throw S("no pointwise operation on curve"); @@ -269,7 +269,7 @@ void NOperate::Integral() if( T->has_dummy() ) throw S("dummy argument t not allowed when operating on data"); - while( ROlo fin=fiter.get() ) { + while( const COlo *fin=fiter.next() ) { int k = fiter.k(); int nz = fin->ZCo.size(); bool savable = nz > 0; @@ -277,10 +277,10 @@ void NOperate::Integral() POld fout; PSpec sout; if ( savable ) { - if ( ROlc fc = R2C(fin) ) - fout = POld( new COld( fc.get() ) ); + if ( const COlc *fc = dynamic_cast<const COlc*>(fin) ) + fout = POld( new COld( fc ) ); else - fout = POld( (R2D(fin))->new_old() ); + fout = POld( (dynamic_cast<const COld*>(fin))->new_old() ); fout->lDoc.push_back( "oi " +expr ); fout->ZCo.pop_back(); fout->xco = fin->ZCo.back(); @@ -349,7 +349,7 @@ void NOperate::IntXY( string mode ) xco = CCoord( "x", "" ); yco = CCoord( "y", "" ); - while( ROld fin = fiter.getD() ) { + while( const COld *fin = fiter.nextD() ) { int nz = fin->ZCo.size(); POld fout( fin->new_old() ); @@ -397,7 +397,7 @@ void NOperate::IntXY( string mode ) void NOperate::Functional( const string& subcmd ) { FileIterator fiter(NSel::selD()); - while ( ROld fin = fiter.getD() ) { + while ( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); fout->lDoc.push_back( "of"+subcmd ); diff --git a/pub/lib/plot.cpp b/pub/lib/plot.cpp index 2c0ae252e043eca4ba51cbfb9e297d1166ce43aa..b114136ab8f5a257c190e135ab4aa5d198e8fec9 100644 --- a/pub/lib/plot.cpp +++ b/pub/lib/plot.cpp @@ -38,7 +38,7 @@ void determine_Xrange( CPlot* plot, vector<int>& JSel ) FileIterator fiter(NSel::selD()); const vector<double> *vx; bool first = true; - while( ROld fd = fiter.getD() ) { + while( const COld *fd = fiter.nextD() ) { // JSel->evaluate( 0, fd->nJ()-1 ); for ( int iv=0; iv<JSel.size(); ++iv ) { vx = &( fd->VS( JSel[iv] )->x ); @@ -68,12 +68,12 @@ void determine_Yrange( CPlot* plot, vector<int>& JSel ) { double inf=INFINITY, sup=-INFINITY; FileIterator fiter(NSel::sel()); - while( ROlo f = fiter.get() ){ + while( const COlo *f = fiter.next() ){ int k = fiter.k(); - ROlc fc = R2C( f ); + const COlc *fc = dynamic_cast<const COlc*>( f ); if ( fc && fc->kconv!=-1 ) continue; - ROld fd = R2D( f ); + const COld *fd = dynamic_cast<const COld*>( f ); // JSel->evaluate( 0, f->nJ()-1 ); for ( int iv=0; iv<JSel.size(); ++iv ){ int j = JSel[iv]; @@ -111,9 +111,9 @@ void determine_Yrange( CPlot* plot, vector<int>& JSel ) //! Plot scan j of data file fd; return number of points plotted. -int plot_data( CPlot* plot, ROld fd, int k, int j, int pstyle ) +int plot_data( CPlot* plot, const COld *fd, int k, int j, int pstyle ) { - const RSpec s = fd->VS(j); + const CSpec *s = fd->VS(j).get(); int n=s->x.size(); if ( n!=s->y.size() ) throw S("BUG: plot: x.size<>y.size"); @@ -179,11 +179,11 @@ int plot_data( CPlot* plot, ROld fd, int k, int j, int pstyle ) //! Plot scan j of convolved curve file fc; return number of points plotted. -int plot_curve_convolved( CPlot* plot, ROlc fc, int k, int j, int cstyle ) +int plot_curve_convolved( CPlot* plot, const COlc *fc, int k, int j, int cstyle ) { vector<double> novec; int kconv, jconv; - RSpec sconv; + const CSpec *sconv; NCurveFile::get_conv( &sconv, &kconv, &jconv, k, j ); vector<double> xp; for ( int i=0; i<sconv->size(); ++i ){ @@ -212,11 +212,11 @@ int plot_curve_convolved( CPlot* plot, ROlc fc, int k, int j, int cstyle ) //! Plot scan j of curve file fc, using grid from file kd. -int plot_curve_to_grid( CPlot* plot, ROlc fc, int k, int j, int cstyle ) +int plot_curve_to_grid( CPlot* plot, const COlc *fc, int k, int j, int cstyle ) { if( fc->kd==-1 ) throw S("data reference not set, incompatible with cg+"); - ROld fd = NOlm::mem_get_D(fc->kd); + const COld *fd = NOlm::mem_get_D(fc->kd); if( !fd ) throw S("data reference not found"); if( j>=fd->nJ() ) @@ -240,7 +240,7 @@ int plot_curve_to_grid( CPlot* plot, ROlc fc, int k, int j, int cstyle ) //! Plot scan j of non-convolved curve file fc on equidistant grid; return number of points plotted. -int plot_curve_equidist( CPlot* plot, ROlc fc, int k, int j, int cstyle ) +int plot_curve_equidist( CPlot* plot, const COlc *fc, int k, int j, int cstyle ) { vector<double> novec; // equidistant grid @@ -264,7 +264,7 @@ int plot_curve_equidist( CPlot* plot, ROlc fc, int k, int j, int cstyle ) //! Plots scan j of non-convolved curve file fc, refining the grid where necessary; return number of points plotted. -int plot_curve_refine( CPlot* plot, ROlc fc, int k, int j, int cstyle ) +int plot_curve_refine( CPlot* plot, const COlc *fc, int k, int j, int cstyle ) { // refinement cares to: // - stop/start line when leaving/reentering/crossing(TODO) the frame @@ -452,10 +452,10 @@ void NPlot::plot( class CPlot *plot, bool add, const string& mode ) } // plot: - while( ROlo f = fiter.get() ){ + while( const COlo *f = fiter.next() ){ int k = fiter.k(); - ROld fd = R2D( f ); - ROlc fc = R2C( f ); + const COld *fd = dynamic_cast<const COld*>( f ); + const COlc *fc = dynamic_cast<const COlc*>( f ); plot->doc_TxLine( f->name ); for (int i=0; i<f->lDoc.size(); i++) diff --git a/pub/lib/ptr.hpp b/pub/lib/ptr.hpp index e23c4d98e5acf1a4fbcaa881ac9ae7dee828369c..ef038ec781b883d0cbc00276a24390ecf0d6bbce 100644 --- a/pub/lib/ptr.hpp +++ b/pub/lib/ptr.hpp @@ -15,12 +15,6 @@ using P = std::shared_ptr< class C>; \ using R = std::shared_ptr<const class C> -DEF_PTR( COlo, POlo, ROlo ); -DEF_PTR( COlc, POlc, ROlc ); -DEF_PTR( COld, POld, ROld ); -DEF_PTR( CSlice, PSlice, RSlice ); -DEF_PTR( CSpec, PSpec, RSpec ); -DEF_PTR( CCurve, PCurve, RCurve ); DEF_PTR( CObj, PObj, RObj ); DEF_PTR( CObjNum, PObjNum, RObjNum ); DEF_PTR( CObjInt, PObjInt, RObjInt ); @@ -36,5 +30,12 @@ DEF_PTR( CObjVecEnu, PObjVecEnu, RObjVecEnu ); #undef DEF_PTR -using RNode = std::shared_ptr<const class CNode>; -using RRef = std::shared_ptr<const class CRef1>; +using POlo = std::shared_ptr<class COlo>; +using POld = std::shared_ptr<class COld>; +using POlc = std::shared_ptr<class COlc>; +using PSlice = std::shared_ptr<class CSlice>; +using PSpec = std::shared_ptr<class CSpec>; +using PCurve = std::shared_ptr<class CCurve>; + +using RNode = std::shared_ptr<const class CNode>; +using RRef = std::shared_ptr<const class CRef1>; diff --git a/pub/lib/reduce_curv.cpp b/pub/lib/reduce_curv.cpp index 0fbb86e7807c912779296fb9ac0d2a538d036620..36a823c4c49e11185bf75522479b3e9be5d4001c 100644 --- a/pub/lib/reduce_curv.cpp +++ b/pub/lib/reduce_curv.cpp @@ -35,7 +35,7 @@ namespace NCvin { double numint_epsabs=1e-10; double numint_epsrel=1e-10; // Auxiliary functions: - double NumericIntegral( ROlc fc, int k, int j, double low, double hig ); + double NumericIntegral( const COlc *fc, int k, int j, double low, double hig ); } @@ -75,7 +75,7 @@ void NCvin::set_integration_tuning_pars( string which ) //! Data transfer between frida and gsl_function. typedef struct { - ROlc fc; + const COlc *fc; int k; int j; } EvalDatTyp; @@ -92,7 +92,7 @@ double myeval( double x, void* data ) //! Numeric integration of one curve. -void cvin_integrate( double *r, double *dr, ROlc fc, int k, int j, const vector<double>& args ) +void cvin_integrate( double *r, double *dr, const COlc *fc, int k, int j, const vector<double>& args ) { EvalDatTyp data; data.fc = fc; diff --git a/pub/lib/reduce_curv.hpp b/pub/lib/reduce_curv.hpp index aae5cfa9f8121f2602df30f379d559791f484af0..f5d4fc973febfa0cb10be87333f490968ee61e98 100644 --- a/pub/lib/reduce_curv.hpp +++ b/pub/lib/reduce_curv.hpp @@ -12,7 +12,7 @@ //! Function type: Evaluate curve functional typedef void (*cvin_eval) ( - double *r, double *dr, ROlc fc, int k, int j, const vector<double>& args ); + double *r, double *dr, const COlc *fc, int k, int j, const vector<double>& args ); //! A wrapper holding a functional. diff --git a/pub/lib/special.cpp b/pub/lib/special.cpp index 48cd9e1844895bec3f8fdc4255f7e90dc10e21c7..334c7ebba45763dc7443b5a98c8b3ca5f72f090e 100644 --- a/pub/lib/special.cpp +++ b/pub/lib/special.cpp @@ -36,7 +36,7 @@ void NSpecial::FourierCosine() int n; FileIterator fiter(NSel::selD()); - while( ROld fin = fiter.getD() ) { + while( const COld *fin = fiter.nextD() ) { POld fout( fin->new_old() ); fout->lDoc.push_back("Fourier Cosine"); diff --git a/pub/lib/toplevel.cpp b/pub/lib/toplevel.cpp index 87187134d9f9d6fda00c2c5929b491f22219503c..673731af1fe0fbc17356a6f2d51b3703aa8cf7b2 100644 --- a/pub/lib/toplevel.cpp +++ b/pub/lib/toplevel.cpp @@ -49,7 +49,7 @@ void fstring_initialize(); // implemented in fstring.cpp; no .hpp file //! Error handler for C code from the GNU Scientific Library. -void my_gsl_error_handler ( const char * reason, const char * file, int line, int gsl_errno ) +void my_gsl_error_handler ( const char * reason, const char *file, int line, int gsl_errno ) { cout<<"GSL error "<<gsl_errno<<": "<<reason<<" in "<<file<<", line "<<line<<"\n"; // throw does not work (probably because GSL is in C, not in C++) @@ -107,7 +107,7 @@ void CFrida::execute_cmd( const string cmdline ) throw S("dummy argument t not allowed outside curve definitions"); if( T->k_dependent() ) { FileIterator fiter(NSel::sel()); - while ( (fiter.get()) ) { + while ( (fiter.next()) ) { int k = fiter.k(); if ( fiter.size()>1 ) cout << "f" << k <<":";