diff --git a/pub/src/expr.cpp b/pub/src/expr.cpp index 3619c7d94eddfe593d725d15fd6e0eee91810c55..6ee18c1ea7ab993d026d86df5d50e97230da618b 100644 --- a/pub/src/expr.cpp +++ b/pub/src/expr.cpp @@ -10,18 +10,15 @@ #include <stdlib.h> // for boost/shared_ptr #include <stdio.h> -#include <ask_simple_value.h> // for calculator -#include <readln.h> +#include <readln.h> // for macro parameters TODO -#include "mystd.h" +#include "mystd.h" // for Integrate #include "olf.h" #include "mem.h" #include "func.h" #include "expr.h" -bool debug = false; -extern bool allow_slow_conv; -#define DEB(args) if ( debug ) printf args +extern bool allow_slow_conv; // TODO unify parameter treatment //***************************************************************************// @@ -776,49 +773,43 @@ void CTree::tree_frf( CResult& ret, const CContext& ctx ) const //* CTree / convolution evaluation *// //***************************************************************************// -//! Retrieve kconv. +//! Return kconc and jconv, according to fc[k,j]->kconv. -uint CTree::setConvK( uint k ) +void CTree::setConv( uint& kconv, uint& jconv, uint k, uint j ) { + // set and check kconv: + if ( k >= NOlm::MOM.size() ) + throw string( "BUG: setConv called with invalid k" ); POlc fc = P2C( NOlm::MOM[k] ); if( !fc ) throw string( "k refers to data, curve expected" ); - return fc->kconv; -} - - -//! Determine jconv. + kconv = fc->kconv; -uint CTree::setConvJ( uint k, uint j, uint kconv ) -{ - POlo f; - POld fconv; - if( kconv==-1 ){ - return -1; - } else { - // check k,j - if ( k >= NOlm::MOM.size() ) - throw string( "invalid k" ); - f = NOlm::MOM[k]; - if ( j>=f->nJ() ) - throw string( "invalid j" ); - // check kconv - if ( kconv >= NOlm::MOM.size() ) - throw string( "invalid kconv" ); - fconv = P2D( NOlm::MOM[kconv] ); - if( !fconv ) - throw string( "cannot convolute with curve" ); - // set jconv - if ( fconv->nJ()==1 ) - return 0; - else if ( fconv->nJ()==f->nJ() ) - return j; - else - throw "cannot convolute file " + strg(k) + " with resolution " + - strg(kconv) + ": number of spectra does not match"; - } + if( kconv==(uint)-1 ){ + jconv = -1; + return; + } else if ( kconv >= NOlm::MOM.size() ) + throw string( "invalid convolution file index " ) + strg(kconv); + + POld fconv = P2D( NOlm::MOM[kconv] ); + if( !fconv ) + throw string( "cannot convolute with curve" ); + + // set and check jconv: + if ( j>=fc->nJ() ) + throw string( "BUG: setConv called with invalid j" ); + + // set jconv: + if ( fconv->nJ()==1 ) // nJ:1 correspondence + jconv = 0; + else if ( fconv->nJ()==fc->nJ() ) // nJ:nJ correspondence + jconv = j; + else + throw "cannot convolute file " + strg(k) + " with resolution " + + strg(kconv) + ": number of spectra does not match"; } + //! Convolution context. Declared locally because only needed in tree_conv. typedef struct { @@ -857,8 +848,8 @@ void CTree::tree_conv( CResult& ret, const CContext& ctx ) const ret.preset_vd( n ); ret.dvd.clear(); // Resolution is given by context 'ctx': - uint kconv = setConvK( ctx.k ); - uint jconv = kconv == -1 ? -1 : setConvJ( ctx.k, ctx.j, kconv ); + uint kconv, jconv; + setConv( kconv, jconv, ctx.k, ctx.j ); double conv_norm; // normalization constant (integral of resolution) double conv_step; // step in x/t, if equidistant PSpec sv; diff --git a/pub/src/expr.h b/pub/src/expr.h index ccdbfb4de761d5f7473bf1aa1e3e7b73e9c17b37..08fb51242782f8695570f89935b0f79a6e4adc62 100644 --- a/pub/src/expr.h +++ b/pub/src/expr.h @@ -140,8 +140,7 @@ class CTree { CTree( char which, const PTree& a0 ); CTree( char which, const PTree& a0, const PTree& a1 ); - static uint setConvK( uint k ); - static uint setConvJ( uint k, uint j, uint kconv ); + static void setConv( uint& kconv, uint& jconv, uint k, uint j ); uint npar() const; diff --git a/pub/src/plot.cpp b/pub/src/plot.cpp index e6a823c85c0d95cdb19d7b7a775b17405c7caae2..aae532a159dad9814aa40c20cfac4ebe4757dce8 100644 --- a/pub/src/plot.cpp +++ b/pub/src/plot.cpp @@ -296,8 +296,8 @@ void NPlot::Plot( class CPlot *plot, bool add ) } } else if ( fc && fc->kconv!=-1 ) { // plot curve, but only at x points of convolutand - uint kconv = fc->kconv; - uint jconv = CTree::setConvJ( k, j, kconv ); + uint kconv, jconv; + CTree::setConv( kconv, jconv, k, j ); POld fconv = P2D( NOlm::MOM[kconv] ); if ( !fconv ) throw string( "convolutand is not a data file" );