Skip to content
Snippets Groups Projects
Commit b496661e authored by Wuttke, Joachim's avatar Wuttke, Joachim
Browse files

correct cauchy; conv debmsg

parent 44233398
No related branches found
No related tags found
No related merge requests found
......@@ -568,6 +568,7 @@ void CTree::tree_val( CTOut *ret, const CContext *ctx ) const
throw string( "2nd arg of conv must be scalar shift" );
double shift = r.d;
if ( kconv==-1 ) {
printf( "DEB no conv %i\n", n );
// Copy convolutand (= unconvoluted theory).
CTOut r;
// unconvoluted theory given by arg[0]:
......@@ -575,6 +576,7 @@ void CTree::tree_val( CTOut *ret, const CContext *ctx ) const
for ( uint i=0; i<n; ++i )
ret->vd[i] = r.to_d(i) - shift;
} else if ( NRead::iofmac("\\cv_mode")==-1 ) {
printf( "DEB conv %i %i step %f\n", n, nv, conv_step );
CTOut r;
CContext cvctx(*ctx);
vector<double> tt;
......@@ -609,8 +611,8 @@ void CTree::tree_val( CTOut *ret, const CContext *ctx ) const
tt[i] = (*(ctx->fargs))[i] - sv->x[iv] - shift;
arg[0]->tree_val( &r, &cvctx );
for ( uint i=0; i<n; ++i )
ret->vd[i] += r.to_d(i) * dx * sv->y[iv] /
conv_norm;
ret->vd[i] += r.to_d(i) * sv->y[iv] *
dx / conv_norm;
}
}
} else {
......@@ -685,6 +687,7 @@ void CTree::tree_val( CTOut *ret, const CContext *ctx ) const
vt = rarg.vd;
}
// evaluate function:
printf( "DEB set FRF context: %i -> %i,%i; x %f .. %f\n", ctx->k, frf.cr.k, frf.cr.j, vt[0], vt[vt.size()-1] );
CContext cu_ctx( frf.cr.k, frf.cr.j );
cu_ctx.set_fargs( &vt );
frf.cr.fc->T->tree_val( ret, &cu_ctx );
......@@ -801,7 +804,6 @@ uint CTree::setConvJ( uint k, uint j, uint kconv ) const
throw "cannot convolute file " + strg(k) + " with resolution " +
strg(kconv) + ": number of scans does not match";
}
}
......
......@@ -117,7 +117,7 @@ double func_gauss (double x, double s) {
double func_gaussnn (double x, double s) { // Gauss not normalized
return s==0 ? 0 : exp(-x*x/(2*s*s)); }
double func_cauchy (double x, double s) {
return s==0 ? 0 : 2*s/twopi/(x*x+s*s); }
return s==0 ? 0 : 2*fabs(s)/twopi/(x*x+s*s); }
double func_kwwc( double v, double a){
if( a<0.1 || a>2.0 ) return 0;
return kwwcf( (float)v, (float)a ); }
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment