fast conv filtering small edits

2d files no longer needed with kiss_fftnd
This commit is contained in:
Mark Borgerding 2004-01-24 01:34:01 +00:00
parent ca97282da6
commit 25e377eaa8
5 changed files with 30 additions and 122 deletions

View File

@ -87,4 +87,4 @@ selftest_short.c:
./mk_test.py -s 10 12 14 > selftest_short.c
clean:
rm -f *~ bm_* st_* tr_* kf_* tkfc_* ff_* ffr_* *.pyc *.pyo
rm -f *~ bm_* st_* tr_* kf_* tkfc_* ff_* ffr_* *.pyc *.pyo *.dat

View File

@ -1,75 +0,0 @@
/*
Copyright (c) 2003, Mark Borgerding
All rights reserved.
Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
* Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
* Neither the author nor the names of any contributors may be used to endorse or promote products derived from this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*/
#include "_kiss_fft_guts.h"
typedef struct {
int minus2; /*magic to signify a 2-d transform*/
kiss_fft_state * rowst;
kiss_fft_state * colst;
kiss_fft_cpx * tmpbuf;
}kiss_fft2d_state;
void * kiss_fft2d_alloc(int nrows,int ncols,int inverse_fft,void*mem,size_t*lenmem)
{
kiss_fft2d_state *st = NULL;
size_t size1, size2, sizetmp, memneeded;
kiss_fft_alloc (ncols, inverse_fft, NULL, &size1);
kiss_fft_alloc (nrows, inverse_fft, NULL, &size2);
sizetmp = sizeof (kiss_fft_cpx) * (ncols > nrows ? ncols : nrows);
memneeded = sizeof (kiss_fft2d_state) + size1 + size2 + sizetmp;
if (lenmem == NULL) {
st = (kiss_fft2d_state *) malloc (memneeded);
} else {
if (*lenmem >= memneeded)
st = (kiss_fft2d_state *) mem;
*lenmem = memneeded;
}
if (!st)
return NULL;
st->minus2 = -2;
st->rowst = (kiss_fft_state *) (st + 1); /*just beyond kiss_fft2d_state struct */
st->colst = (kiss_fft_state *) ((char *) (st->rowst) + size1);
st->tmpbuf = (kiss_fft_cpx *) ((char *) (st->rowst) + size1 + size2);
kiss_fft_alloc (ncols, inverse_fft, st->rowst, &size1);
kiss_fft_alloc (nrows, inverse_fft, st->colst, &size2);
return st;
}
void kiss_fft2d(const void * cfg,const kiss_fft_cpx *fin,kiss_fft_cpx *fout)
{
/* input buffer fin is stored row-wise */
kiss_fft2d_state *st = ( kiss_fft2d_state *)cfg;
int row,col;
int nrows,ncols;
nrows = st->colst->nfft;
ncols = st->rowst->nfft;
/*fft each column*/
for (col=0;col<ncols;++col) {
for (row=0;row< nrows ;++row)
st->tmpbuf[row] = fin[row*ncols + col];
kiss_fft(st->colst,st->tmpbuf,st->tmpbuf);
for (row=0;row< nrows ;++row) {
fout[row*ncols + col] = st->tmpbuf[row];
}
}
/*fft each row */
for (row=0;row< nrows ;++row)
kiss_fft(st->rowst , fout + row*ncols , fout + row*ncols );
}

View File

@ -1,12 +0,0 @@
#ifndef KISS_FFT2D_H
#define KISS_FFT2D_H
#include "kiss_fft.h"
/* allocate a 2-dimensional FFT
the data should be stored rowwise,
in other words, an array made up of row[0], then row[1], etc
*/
void * kiss_fft2d_alloc(int nrows,int ncols,int inverse_fft,void * mem,size_t * lenmem);
void kiss_fft2d(const void* cfg_from_alloc , const kiss_fft_cpx *fin,kiss_fft_cpx *fout );
#endif

View File

@ -87,4 +87,4 @@ selftest_short.c:
./mk_test.py -s 10 12 14 > selftest_short.c
clean:
rm -f *~ bm_* st_* tr_* kf_* tkfc_* ff_* ffr_* *.pyc *.pyo
rm -f *~ bm_* st_* tr_* kf_* tkfc_* ff_* ffr_* *.pyc *.pyo *.dat

View File

@ -30,7 +30,7 @@ typedef kiss_fft_cpx kffsamp_t;
static int verbose=0;
void * kiss_fastfir_alloc(const kffsamp_t * imp_resp,size_t n_imp_resp,
int nfft,void * mem,size_t*lenmem);
size_t * nfft,void * mem,size_t*lenmem);
typedef struct {
int nfft;
@ -47,7 +47,7 @@ static const kiss_fft_cpx CZERO={0,0};
void * kiss_fastfir_alloc(
const kffsamp_t * imp_resp,size_t n_imp_resp,
int nfft, /* if <= 0, an appropriate size will be chosen */
size_t *pnfft, /* if <= 0, an appropriate size will be chosen */
void * mem,size_t*lenmem)
{
kiss_fastfir_state *st = NULL;
@ -55,8 +55,11 @@ void * kiss_fastfir_alloc(
size_t memneeded = sizeof(kiss_fastfir_state);
char * ptr;
size_t i;
int nfft=0;
float scale;
int n_freq_bins;
if (pnfft)
nfft=*pnfft;
if (nfft<=0) {
/* determine fft size as next power of two at least 2x
@ -79,10 +82,8 @@ void * kiss_fastfir_alloc(
/*ifftcfg*/
FFT_ALLOC (nfft, 1, NULL, &len_ifftcfg);
memneeded += len_ifftcfg;
/* tmpbuf */
memneeded += sizeof(kffsamp_t) * nfft;
/* fir_freq_resp */
memneeded += sizeof(kiss_fft_cpx) * n_freq_bins;
/* freqbuf */
@ -123,7 +124,7 @@ void * kiss_fastfir_alloc(
memset(st->tmpbuf,0,sizeof(kffsamp_t)*nfft);
/*zero pad in the middle to left-rotate the impulse response
* This puts the scrap samples at the end of the inverse fft'd buffer */
This puts the scrap samples at the end of the inverse fft'd buffer */
st->tmpbuf[0] = imp_resp[ n_imp_resp - 1 ];
for (i=0;i<n_imp_resp - 1; ++i) {
st->tmpbuf[ nfft - n_imp_resp + 1 + i ] = imp_resp[ i ];
@ -138,6 +139,8 @@ void * kiss_fastfir_alloc(
st->fir_freq_resp[i].r *= scale;
st->fir_freq_resp[i].i *= scale;
}
if (pnfft)
*pnfft = nfft;
return st;
}
@ -145,7 +148,6 @@ static void fastconv1buf(const kiss_fastfir_state *st,const kffsamp_t * in,kffsa
{
int i;
FFTFWD( st->fftcfg, in , st->freqbuf );
/* multiply the frequency response of the input signal by*/
/* that of the fir filter*/
for ( i=0; i<st->n_freq_bins; ++i ) {
@ -182,8 +184,6 @@ size_t kff_nocopy(
#ifdef FAST_FILT_UTIL
#define BUFLEN 1024
void do_filter(
FILE * fin,
FILE * fout,
@ -191,53 +191,48 @@ void do_filter(
size_t n_imp_resp,
size_t nfft)
{
void * cfg = kiss_fastfir_alloc(imp_resp,n_imp_resp,nfft,0,0);
size_t max_inbuf=5*nfft;
size_t max_outbuf=5*nfft;
kffsamp_t *inbuf;
kffsamp_t *outbuf;
void * cfg = kiss_fastfir_alloc(imp_resp,n_imp_resp,&nfft,0,0);
size_t max_buf=5*nfft;
kffsamp_t *inbuf = (kffsamp_t*)malloc(max_buf*sizeof(kffsamp_t));
kffsamp_t *outbuf = (kffsamp_t*)malloc(max_buf*sizeof(kffsamp_t));
size_t ninbuf,noutbuf;
size_t idx_inbuf=0;
int done=0;
inbuf = (kffsamp_t*)malloc(max_inbuf*sizeof(kffsamp_t));
outbuf = (kffsamp_t*)malloc(max_outbuf*sizeof(kffsamp_t));
do{
int zpad=0;
ninbuf = fread(&inbuf[idx_inbuf],sizeof(inbuf[0]),max_inbuf-idx_inbuf,fin );
do{
ninbuf = fread(&inbuf[idx_inbuf],sizeof(inbuf[0]),max_buf-idx_inbuf,fin );
if (ninbuf==0) {
/* zero pad the input to the fft size */
done=1;
zpad = nfft - idx_inbuf;
memset(&inbuf[idx_inbuf],0,sizeof(inbuf[0])*zpad);
if (verbose) fprintf(stderr,"zero padding %d samples\n",zpad);
ninbuf = nfft;
/* We don't strictly need to do the zero padding,since the output samples
* will get thrown away anyway. But the zeros should help keep the fft magnitudes in
* check and thus decrease errors. */
memset(&inbuf[idx_inbuf],0,sizeof(inbuf[0])*zpad);
ninbuf = nfft;/* force an fft */
}else{
ninbuf += idx_inbuf;
}
/*perform the fast convolution filtering*/
noutbuf = kff_nocopy( cfg, inbuf, outbuf, ninbuf );
if (verbose) fprintf(stderr,"kff_nocopy(,,,%d) -> %d\n",ninbuf,noutbuf);
/* move the unconsumed samples to the front */
idx_inbuf = ninbuf-noutbuf;
/* move the unconsumed input samples to the front */
idx_inbuf = ninbuf - noutbuf;
memmove( inbuf , &inbuf[noutbuf] , sizeof(inbuf[0])*(ninbuf-noutbuf) );
if (verbose) fprintf(stderr,"moved %d samples to front of buffer\n",idx_inbuf);
/*whatever padding was done to get the size up to nfft, must now be ignored*/
noutbuf -= zpad;
if ( fwrite( outbuf,
sizeof(outbuf[0]),
noutbuf,
fout)
!= noutbuf )
{
if ( fwrite( outbuf, sizeof(outbuf[0]), noutbuf, fout) != noutbuf ) {
fprintf(stderr,"short write %d \n",noutbuf);
fprintf(stderr,"zpad= %d \n",zpad);
exit(1);
}
}while ( ! done );
}while ( ! zpad );
fclose(fout);
free(cfg);
free(inbuf);