estimation.h

软件无线电的平台

/**
 * Channel Estimation Routines
 * 
 * Alessandro Nordio, 01/03
 * Linus Gasser 02/10/25
 *
 * This version is optimized for MMX.
 *
 */

#include "mmx.h"


/**
 * @author Alessandro Nordio
 *
 * 23/03/2001
 *
 *  
 * @param x input vector of 2048 short elements, is assumed to be real 
 *  and in the format |x0 x1 x2 x3 x4 ... x(N-1) 
 * @param seq is the pointer to the FFT of the midamble sequence 
 *  (already precomputed and in the format)
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'seq' contains 512 mmx_t elements 
 * @param tw twiddle factor for FFT in the format 
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'tw' contains 2047 elements 
 *  ordered as described in the function fft_radix2_fixpoint_mmx
 * @param tw2 twiddle factor for IFFT in the format 
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'tw2' contains 2047 elements 
 *  ordered as described in the function fft_radix2_fixpoint_mmx
 * @param rev pointer to the reverse buffer. It contains 2048 
 *  short unsigned integers
 */
int channel_estimation_2048(mmx_t *x, mmx_t *y, mmx_t *seq,  
			    mmx_t *tw, mmx_t *tw2, unsigned short *rev);


/**
 * @author Alessandro Nordio
 *
 * 23/03/2001
 *
 *  
 * @param x input vector of 512 short elements, is assumed to be real 
 *  and in the format |x0 x1 x2 x3 x4 ... x(N-1) 
 * @param seq is the pointer to the FFT of the midamble sequence 
 *  (already precomputed and in the format)
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'seq' contains 128 mmx_t elements 
 * @param tw twiddle factor for FFT in the format 
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'tw' contains 511 elements 
 *  ordered as described in the function fft_radix2_fixpoint_mmx
 * @param tw2 twiddle factor for IFFT in the format 
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'tw2' contains 511 elements 
 *  ordered as described in the function fft_radix2_fixpoint_mmx
 * @param rev pointer to the reverse buffer. It contains 512
 *  short unsigned integers
 */
int channel_estimation_512(mmx_t *x,mmx_t *y, mmx_t *seq,  mmx_t *tw, 
			   mmx_t *tw2, unsigned short *rev);

/**
 * @author Alessandro Nordio
 *
 * 23/03/2001
 * Channel Estimation routine for 3GPP midamble type 2 and 
 * oversampling 4 (SAMPLES PER CHIP=4) (Midable length [chips] = 192
 * this function uses a 192*4 = 256*3 = 768 point FFT and IFFT )
 *  
 * @param x input vector of 768 short elements, is assumed to be real 
 *  and in the format |x0 x1 x2 x3 x4 ... x(N-1) 
 * @param seq is the pointer to the FFT of the midamble sequence 
 *  (already precomputed and in the format)
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'seq' contains 192 mmx_t elements 
 * @param tw twiddle factor for FFT in the format 
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'tw' contains 1791 elements 
 *  ordered as described in the function fft_radix2_fixpoint_mmx
 * @param tw2 twiddle factor for IFFT in the format 
 *  |Re0 -Im0 Im0 Re0 |Re1 -Im1 Im1 Re1 | ... 'tw2' contains 1791 elements 
 *  ordered as described in the function fft_radix2_fixpoint_mmx
 * @param rev pointer to the reverse buffer. It contains 768
 *  short unsigned integers given by the function reverse_buffer768
 */
int channel_estimation_768(mmx_t *x,mmx_t *y, mmx_t *seq, mmx_t *tw, 
			   mmx_t *tw2, unsigned short *rev);


/**
 * This function computes a 768 point complex FFT in fixed point
 *
 * @author Alessandro Nordio 22/03/2001
 *
 * @param x input vector of 768 mmx_t elements with format 
 *  |Re0,Im0,Re0,Im0|......|Re(N-1),Im(N-1),Re(N-1),Im(N-1)|
 * @param y output vector with N elements and the same format 
 *  of the input. The input vector can be reused for the output. 
 *  (e.g.  x = &t[0], y = &t[0]; x and y points to the same vector)
 * @param tw the twiddle here is a vector of 255+2*768 mmx_t 
 *  elements. Is a concatenation of 3 parts:<ul>
 *  <li>the first  part of 255 elements contains the twiddles for a 
 *   256 point fft organized as described in the function 
 *   'fft_radix2_fixpoint_mmx' (see comments insife the function 
 *   and consider log2size = 8)</li>
 *  <li>the second part of 768 elements contains exp(-j*2*pi*n/N)   
 *   for n=0..(N-1) in the format |Re,-Im,Im,Re|</li>
 *  <li>the third  part of 768 elements contains exp(-j*2*pi*2*n/N) 
 *   for n=0..(N-1) in the format |Re,-Im,Im,Re|</li>
 * </ul>
 * @param rev pointer to the reverse buffer
 * @param log2_amp this is an unsigned integer that for default is 0.
 *  when not 0 it increases the output dinamic by a factor 2^log2_amp
 * ** WARNING: values different from 0 can cause overflow in the FFT!!!! **
 * ** WARNING: log2_amp must be lower than 7
 * @param flag<ul>
 *  <li>0 the input is real, computes only the first half of the 
 *   output vector (see butterfly3)</li>
 *  <li>1 the input is in mmx_t complex format</li>
 *  <li>2 the input is as in 1 but the second half of x is assumed 
 *   to be 0 (avoids first IFFT stage) computes the overall output vector</li>
 * </ul>
 */
int fft_768_fixpoint_mmx(mmx_t *x, mmx_t *y, mmx_t *twiddle, 
			 unsigned short *rev, unsigned short log2_amp,
			 unsigned short flag);


/*
 * flag = 0 input is real x0,x1,x2,x3,...
 * flag = 1 input is cpx and in mmx_t format Re0 Im0 re0 Im0, ....
 * flag = 2 as in 1 but the second half of the input vector is assumed to be 0
 */
int reverse_order768(mmx_t *x,mmx_t *y,unsigned short *rev,unsigned short flag);

/**
 * compute the reverse budder for a 768 point FFT
 * @param buf must point to a 768 short (16 bit) element vector
 */
int reverse_buffer768(unsigned short *buf);


/**
 * @short This function computes N butterflies of size 3
 *
 * @author Alessandro Nordio 22/03/2001
 *
 * @param x input vector of N mmx_t elements with format 
 *  |Re0,Im0,Re0,Im0| ... |Re(N-1),Im(N-1),Re(N-1),Im(N-1)|
 * @param y output vector with N elements and the same format of the 
 *  input. The input vector can be reused for the output. 
 *  (e.g.  x = &t[0], y = &t[0]; x and y points to the same vector)
 * @param tw the twiddle here is a vector of 2*N mmx_t elements. 
 *  Is a concatenation of two equal length parts:<ul>
 *   <li>the first  part of N elements contains exp(-j*2*pi*n/N)   
 *       for n=0..(N-1) in the format |Re,-Im,Im,Re|</li>
 *   <li>the second part of N elements contains exp(-j*2*pi*2*n/N) 
 *       for n=0..(N-1) in the format |Re,-Im,Im,Re|</li>
 * @param N the size of the input vector (and also the number of 
 *  size-3 butterflies computed)
 * @param flag<ul>
 *  <li>0 computes N regular butterflies</li>
 *  <li>1 computes the output of the butterflies corresponding to the 
 *   first half of the output vector</li>
 *  </ul>
 */
int butterfly3( mmx_t *x, mmx_t *y, mmx_t *tw, 
		unsigned int N, unsigned short flag );


/**
 * @short FFT base 2 in fixed point
 * @author Alessandro Nordio, 03/2001
 *
 * This version is optimized for MMX.  (Version 4.0)
 *
 * @param x input vector (the input vector format depends on 'flag' 
 *  see 'flag' here below)
 * @param y output vector of size N = 2^log2size mmx_t element packed 
 *  in the format [Re0,Im0,Re0,Im0, Re1,Im1,Re1,Im1, ..., 
 *  Re(N-1),Im(N-1),Re(N-1),Im(N-1)]
 * @param twiddle twiddle factors are precomputed and are in the format 
 *  Twiddle format [Re0,-Im0,Im0,-Re0|Re1,-Im1,Im1,-Re1|,.....]
 *  Each element of twiddle is a mmx_t union (64 bit) that contains |Re,-Im,Im,Re|
 *  of a certain twidde factor. 
 *  This way of storing the twiddles is useful to increase the performance of the FFT
 *  The twiddle vector is logically divided into 'log2size' consecutive parts: 
 *  (e.g. log2size=11 for a 2048 point fft).
 *  The i-th part contains n_i=2^i elements, for i=0..log2size-1 .
 *  The elements of the i-th part are  round(A*exp(-j*2*pi*n/(2^(i+1)))) for n=0..n_i
 *  where A is a power-of-two constant integer. We assume A=2^14.
 *  (Since A=2^14 then in the butterflies we shift the output of the mpy by 14 
 *   (see below))
 *  The elements of the twiddle are precomputed using a matlab file 
 *   'twiddle_mmx.m' and are stored in a binary file in little endian format.
 *  The overall twiddle vector contains 2^log2size -1 complex elements.
 *  The i-th part of the twiddle vector is used during the i-th stage of the FFT
 *  In this version the first stage of the FFT is simplified and has been 
 *  taken out of the loop.
 * @param rev reverse buffer. Contains the reverse indexes (N indexes) 
 * @param log2size is the logarithm in base 2 of the size of the FFT.
 *  WARNING log2size>=2
 * @param log2_amp this is an unsigned integer that for default is 0.
 * when not 0 it increases the output dinamic by a factor 2^log2_amp
 * WARNING: values different from 0 can cause overflow in the FFT!!!!
 * WARNING: log2_amp must be lower than log2size-1
 * @param flag flag specifies the input format of the vector x<ul>
 *  <li>0 x contains N = 2^log2size short integers elements (real) 
 *   (Input format:  [x_0,x_1,x_2,x_3,.....x_N] where x_i is a 16 bit integer)
 *   before doing FFT the function does reverse bit ordering</li>
 *  <li>1 x contains N = 2^log2size mmx_t elements 
 *   (Input format: [Re0,Im0,Re0,Im0, Re1,Im1,Re1,Im1 , ... , 
 *   Re(N-1),Im(N-1),Re(N-1),Im(N-1)] before doing FFT the function does 
 *   reverse bit ordering</li>
 *  <li>2 as for 0 but does not do reverse bit ordering (NOT IMPLEMENTED !!)</li>
 *  <li>3 as for 1 but does not do reverse bit ordering</li>
 *  <li>4 as for 1 but the second half of the input vector is assumed 
 *   to be 0 (useful for IFFT in 3GPP channel estimation)
 *   IMPORTANT: the first FFT stage is integrated in the 
 *   reverse bit ordering loop</li>
 *  <li>5 does not do reverse bit ordering and start from the second FFT stage
 *   IMPORTANT: useful for a 768 point IFFT where reordering has been 
 *   done outside the function</li>
 *  </ul>
 */
int fft_radix2_fixpoint_mmx(mmx_t *x, 
			    mmx_t *y,
			    mmx_t *twiddle,
			    unsigned short *rev, 
			    unsigned short log2size, 
			    unsigned short log2_amp, 
			    short flag);


unsigned short ReverseBits ( unsigned index, unsigned NumBits );

int mult_cpx_vector(mmx_t *x1, mmx_t *x2, mmx_t *y, 
		    unsigned int N, unsigned short log2_amp);