diff options
Diffstat (limited to 'editor/specialty/pgmabel.c')
| -rw-r--r-- | editor/specialty/pgmabel.c | 13 |
1 files changed, 6 insertions, 7 deletions
diff --git a/editor/specialty/pgmabel.c b/editor/specialty/pgmabel.c index aa748f88..3e05be5e 100644 --- a/editor/specialty/pgmabel.c +++ b/editor/specialty/pgmabel.c @@ -100,25 +100,24 @@ abel ( float *y, int N, double *adl) { register int n; double *rho, *rhop; /* results and new index */ - float *yp; /* new indizes for the y-array */ MALLOCARRAY(rho, N); if( !rho ) pm_error( "out of memory" ); + for (n=0 ; n<N ; n++) + rho[n] = 0; + rhop = rho; - yp = y; for (n=0 ; n<N ; n++) { - *(rhop++) = ((*yp++) - Sum(n,rho,N,adl))/(adl[n*N+n]); -/* *(rhop++) = ((*yp++) - Sum(n,rho,N))/(dr(n,n+0.5,N)); old version */ - if ( *rhop < 0.0 ) *rhop = 0.0; /* error correction ! */ -/* if (n > 2) rhop[n-1] = (rho[n-2]+rho[n-1]+rho[n])/3.0; stabilization*/ + rhop[n] = MAX(0, (y[n] - Sum(n,rho,N,adl))/(adl[n*N+n])); + /* Clip to 0 for error correction ! */ } for (n=0 ; n<N ; n++) { if (( n>=1 )&&( n<N-1 )) - (*y++) = ((rho[n-1]*0.5+rho[n]+rho[n+1]*0.5)/2.0);/*1D median filter*/ + (*y++) = ((rho[n-1]*0.5+rho[n]+rho[n+1]*0.5)/2.0);/*1D median filter*/ else (*y++) = rho[n]; } free(rho); |