/* -*- C -*- */
/*
 * Copyright (c) 2003 Jindrich Novy (jnovy@users.sf.net)
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 *
 */

void X(convolve_1d_adapt)( FCE_SAMPLE *a1, FCE_SAMPLE *a2, int xs, FCE_FLOAT *v, int o, FCE_FLOAT *l, int ls ) {
	int x, xx, vp, o2 = o>>1;
	FCE_FLOAT s, n, g;
	
	ls>>=1;
	for ( x=0; x<o2; x++ ) {
		for ( n=s=0., vp=o2-x, xx=0; vp<o; xx++, vp++  ) {
			g = v[vp]*l[(int)(ls+a1[x]-a1[xx])];
			s += g*a1[xx];
			n += g;
		}
		a2[x] = s/n;
	}

	for ( ; x<xs-o2; x++ ) {
		for ( n=s=0., vp=0, xx=x-o2; vp<o; xx++, vp++ ) {
			g = v[vp]*l[(int)(ls+a1[x]-a1[xx])];
			s += g*a1[xx];
			n += g;
		}
		a2[x] = s/n;
	}
	
	for ( ; x<xs; x++ ) {
		for ( n=s=0., vp=0, xx=x-o2; xx<xs; xx++, vp++ ) {
			g = v[vp]*l[(int)(ls+a1[x]-a1[xx])];
			s += g*a1[xx];
			n += g;
		}
		a2[x] = s/n;
	}
}

void X(convolve_1d_adapt_stride)( FCE_SAMPLE *a1, FCE_SAMPLE *a2, int xs, FCE_FLOAT *v, int o, FCE_FLOAT *l, int ls, int stride ) {
	int x, xx, vp, o2 = o>>1, x_stride;
	FCE_FLOAT s, n, g;

	ls>>=1;
	for ( x_stride=x=0; x<o2; x++, x_stride+=stride ) {
		for ( n=s=0., vp=o2-x, xx=0; vp<o; xx+=stride, vp++  ) {
			g = v[vp]*l[(int)(ls+a1[x_stride]-a1[xx])];
			s += g*a1[xx];
			n += g;
		}
		a2[x_stride] = s/n;
	}

	for ( n=0., vp=0; vp<o; vp++ ) n += v[vp];

	for ( ; x<xs-o2; x++, x_stride+=stride ) {
		for ( n=s=0., vp=0, xx=(x-o2)*stride; vp<o; xx+=stride, vp++ ) {
			g = v[vp]*l[(int)(ls+a1[x_stride]-a1[xx])];
			s += g*a1[xx];
			n += g;
		}
		a2[x_stride] = s/n;
	}
	
	for ( ; x<xs; x++, x_stride+=stride ) {
		for ( n=s=0., vp=0, xx=(x-o2)*stride; xx<xs*stride; xx+=stride, vp++ ) {
			g = v[vp]*l[(int)(ls+a1[x_stride]-a1[xx])];
			s += g*a1[xx];
			n += g;
		}
		a2[x_stride] = s/n;
	}
}

void X(convolve_2d_adapt)( FCE_SAMPLE *a1, FCE_SAMPLE *a2, FCE_SAMPLE *at, int xs, int ys, FCE_FLOAT *v, int o, FCE_FLOAT *l, int ls ) {
	int x, y, p, at_allocated = 0;

	if ( !at ) { 
		assert( ((void*)at=malloc(xs*ys*sizeof(FCE_SAMPLE))) );
		at_allocated = 1;
	}
	
	for ( p=y=0; y<ys; y++, p+=xs )
		X(convolve_1d_adapt)( &a1[p], &at[p], xs, v, o, l, ls );

	for ( x=0; x<xs; x++ )
		X(convolve_1d_adapt_stride)( &at[x], &a2[x], ys, v, o, l, ls, xs );
		
	if ( at_allocated ) free( at );
}

void X(convolve_3d_adapt)( FCE_SAMPLE *a1, FCE_SAMPLE *a2, FCE_SAMPLE *at, int xs, int ys, int zs, FCE_FLOAT *v, int o, FCE_FLOAT *l, int ls ) {
	int x, y, z, py, pz, at_allocated = 0;

	if ( !at ) { 
		assert( ((void*)at=malloc(xs*ys*zs*sizeof(FCE_SAMPLE))) );
		at_allocated = 1;
	}
	
	for ( pz=z=0; z<zs; z++, pz+=xs*ys ) {
		for ( py=y=0; y<ys; y++, py+=xs )
			X(convolve_1d_adapt)( &a1[pz+py], &a2[pz+py], xs, v, o, l, ls );

		for ( x=0; x<xs; x++ )
			X(convolve_1d_adapt_stride)( &a2[pz+x], &at[pz+x], ys, v, o, l, ls, xs );
	}

	if ( z>1 ) {
		for ( z=0; z<xs*ys; z++ )
			X(convolve_1d_adapt_stride)( &at[z], &a2[z], zs, v, o, l, ls, xs*ys );
	}
	
	if ( at_allocated ) free( at );
}