?? flat_rayleigh.cpp
字號:
/*****************************************************************************
* Author: Christos Komninakis *
* Date: January 4, 2000. *
* Content: C++ code file for the CLASS "flat_rayleigh". *
*****************************************************************************/
#include "flat_rayleigh.h"
#define PI 3.14159265358979
/**************** Constructor ***************************************/
flat_rayleigh::flat_rayleigh(int seeed, double fD, double pwr, bool flag_indep)
{
chan_seed = seeed; //-23;
PWR = pwr;
IndepFlag = flag_indep;
chan_val = Complex(0.0, 0.0);
K = 7; H = 7; H2 = 2*H;
if (fD > 0.2) {
cout << "Warning: Discrete Doppler fDT > 0.2, handled as fDT=0.2 exactly" << endl;
I = 1;
} else
I = (int) (0.2 / fD);
last_i = 0; IP = 0;
int i, t, j, k;
a = new double [K];
b = new double [K];
c = new double [K];
d = new double [K];
double A_o;
st = new Complex* [3];
for (i=0; i<3; i++)
st[i] = new Complex [K];
sinc_matrix = new double* [I];
for (i=0; i<I; i++)
sinc_matrix[i] = new double [H2];
buff_f = new Complex [H2];
/* Following are the values of a designed ROOT DOPPLER filter prototype,
consisting of K=7 biquads, with fD = 0.2 */
a[0] = 0.98323218491754; b[0] = 0.13953807910136;
c[0] = -0.60134822489663; d[0] = 0.48938674026073;
a[1] = -0.58056307877395; b[1] = 0.99977164701498;
c[1] = -0.42846077703117; d[1] = 0.24612477470991;
a[2] = -0.50498097209321; b[2] = 0.99986531836969;
c[2] = -0.63050157648166; d[2] = 0.9350218659307;
a[3] = 0.20534422769486; b[3] = 0.98691526068222;
c[3] = -0.096351029308223; d[3] = -0.060828833103004;
a[4] = 1.2667005604552; b[4] = 0.93699896920708;
c[4] = -0.64949215100022; d[4] = 0.80427856161307;
a[5] = -0.60545728027338; b[5] = 0.99998818741754;
c[5] = -0.61939958040555; d[5] = 0.99895172368533;
a[6] = -0.29077222501287; b[6] = 0.99832970528904;
c[6] = -0.62124137243814; d[6] = 0.98111826315356;
A_o = 0.028241759549181;
/* now really start */
Ao = PWR * A_o / sqrt(2.0);
/* produce the table storing all the values of the interpolating function,
windowed and in reversed order, so that filtering is straightforward */
double scale = 0.0;
for (i=0; i<I; i++) {
for (j=0; j<H2; j++) {
if ((i == 0) & (j == H))
sinc_matrix[i][j] = 1.0;
else {
sinc_matrix[i][j] =
sin(PI*(j-H)+(PI/I)*i) / (PI*(j-H)+(PI/I)*i) * //;
(0.5 - 0.5*cos(2.0*PI*(i+I*j)/(H2*I))); // windowed, Hanning
//(0.54 - 0.46*cos(2.0*PI*(i+I*j)/(H2*I))); // windowed, Hamming
//(0.42 - 0.5*cos(2.0*PI*(i+I*j)/(H2*I)) + 0.08*cos(4.0*PI*(i+I*j)/(H2*I))); // windowed, Blackman
}
scale += (sinc_matrix[i][j])*(sinc_matrix[i][j]);
}
}
scale = sqrt(1.0/scale*I);
// scale and swap
double temp;
for (i=0; i<I; i++) {
for (j=0; j<H; j++) {
temp = sinc_matrix[i][j];
sinc_matrix[i][j] = scale * sinc_matrix[i][H2-1-j];
sinc_matrix[i][H2-1-j] = scale * temp;
}
}
sinc_matrix[0][H] *= scale;
/* work the channel for some time, to overcome transients */
for (t=0; t<I*4*(4*K+H2); t++) {
if (last_i == I-1) {
last_i = 0;
/* now produce new element, and insert it into the buffer: buff_f[IP] */
st[0][0] = Complex(Gaussian(&chan_seed), Gaussian(&chan_seed))*Ao -
st[1][0]*c[0] - st[2][0]*d[0];
for (k=1; k<K; k++)
st[0][k] = st[0][k-1] + st[1][k-1]*a[k-1] + st[2][k-1]*b[k-1] -
st[1][k]*c[k] - st[2][k]*d[k];
buff_f[IP] = st[0][K-1] + st[1][K-1]*a[K-1] + st[2][K-1]*b[K-1];
for (k=0; k<K; k++) {
st[2][k] = st[1][k];
st[1][k] = st[0][k];
}
/* done, produced buff_f[IP], now move IP to next position */
IP = (IP + 1) % H2;
} else
last_i++;
// stop here, no need to interpolate, since idle run
}
}
/******* Function pass_through, without providing any CSI ************/
void flat_rayleigh::pass_through(int length, Complex *inp, Complex *outp)
{
register int k, j, t;
if (IndepFlag) {
chan_seed += 1; // reset the seed
// clean-up and then run idle, to ensure transients are dead
for (k=0; k<K; k++) {
st[0][k] = Complex(0.0, 0.0);
st[1][k] = Complex(0.0, 0.0);
st[2][k] = Complex(0.0, 0.0);
}
for (j=0; j<H2; j++)
buff_f[j] = Complex(0.0, 0.0);
// and now run idle for a little...
for (t=0; t<I*4*(4*K+H2); t++) {
if (last_i == I-1) {
last_i = 0;
/* produce new element, and insert it into buff_f[IP] */
st[0][0] = Complex( Gaussian(&chan_seed), Gaussian(&chan_seed) );
st[0][0] = st[0][0]*Ao - st[1][0]*c[0] - st[2][0]*d[0];
for (k=1; k<K; k++)
st[0][k] = st[0][k-1] + st[1][k-1]*a[k-1] + st[2][k-1]*b[k-1] -
st[1][k]*c[k] - st[2][k]*d[k];
buff_f[IP] = st[0][K-1] + st[1][K-1]*a[K-1] + st[2][K-1]*b[K-1];
for (k=0; k<K; k++) {
st[2][k] = st[1][k];
st[1][k] = st[0][k];
}
/* done, produced buff_f[IP], now move IP to next position */
IP = (IP + 1) % H2;
} else
last_i++;
// end here, don't produce chan_val, since idle run
}
} /* end of dry run for independence */
/* run for real now */
for (t=0; t<length; t++) {
if (last_i == I-1) {
last_i = 0;
/* now produce new element, and insert it into the buffer: buff_f[IP] */
st[0][0] = st[0][0] = Complex( Gaussian(&chan_seed), Gaussian(&chan_seed) );
st[0][0] = st[0][0]*Ao - st[1][0]*c[0] - st[2][0]*d[0];
for (k=1; k<K; k++)
st[0][k] = st[0][k-1] + st[1][k-1]*a[k-1] + st[2][k-1]*b[k-1] -
st[1][k]*c[k] - st[2][k]*d[k];
buff_f[IP] = st[0][K-1] + st[1][K-1]*a[K-1] + st[2][K-1]*b[K-1];
for (k=0; k<K; k++) {
st[2][k] = st[1][k];
st[1][k] = st[0][k];
}
/* done, produced buff_f[IP], now move IP to next position */
IP = (IP + 1) % H2;
} else
last_i++;
buff_sinc = sinc_matrix[last_i];
chan_val = Complex(0.0, 0.0);
for (j=0; j<H2; j++)
chan_val += buff_f[(IP + j) % H2] * buff_sinc[j];
outp[t] = inp[t] * chan_val;
}
}
/******* Overloaded Function pass_through, providing CSI also ************/
void flat_rayleigh::pass_through(int length, Complex *inp,
Complex *outp, Complex *csi)
{
register int k, j, t;
if (IndepFlag) {
chan_seed += 1; // reset the seed
// clean-up and then run idle, to ensure transients are dead
for (k=0; k<K; k++) {
st[0][k] = Complex(0.0, 0.0);
st[1][k] = Complex(0.0, 0.0);
st[2][k] = Complex(0.0, 0.0);
}
for (j=0; j<H2; j++)
buff_f[j] = Complex(0.0, 0.0);
// and now run idle for a little...
for (t=0; t<I*4*(4*K+H2); t++) {
if (last_i == I-1) {
last_i = 0;
/* produce new element, and insert it into buff_f[IP] */
st[0][0] = Complex( Gaussian(&chan_seed), Gaussian(&chan_seed) );
st[0][0] = st[0][0]*Ao - st[1][0]*c[0] - st[2][0]*d[0];
for (k=1; k<K; k++)
st[0][k] = st[0][k-1] + st[1][k-1]*a[k-1] + st[2][k-1]*b[k-1] -
st[1][k]*c[k] - st[2][k]*d[k];
buff_f[IP] = st[0][K-1] + st[1][K-1]*a[K-1] + st[2][K-1]*b[K-1];
for (k=0; k<K; k++) {
st[2][k] = st[1][k];
st[1][k] = st[0][k];
}
/* done, produced buff_f[IP], now move IP to next position */
IP = (IP + 1) % H2;
} else
last_i++;
// end here, don't produce chan_val, since idle run
}
} /* end of dry run for independence */
for (t=0; t<length; t++) {
if (last_i == I-1) {
last_i = 0;
/* now produce new element, and insert it into the buffer: buff_f[IP] */
st[0][0] = Complex( Gaussian(&chan_seed), Gaussian(&chan_seed) )*Ao -
st[1][0]*c[0] - st[2][0]*d[0];
for (k=1; k<K; k++)
st[0][k] = st[0][k-1] + st[1][k-1]*a[k-1] + st[2][k-1]*b[k-1] -
st[1][k]*c[k] - st[2][k]*d[k];
buff_f[IP] = st[0][K-1] + st[1][K-1]*a[K-1] + st[2][K-1]*b[K-1];
for (k=0; k<K; k++) {
st[2][k] = st[1][k];
st[1][k] = st[0][k];
}
/* done, produced buff_f[IP], now move IP to next position */
IP = (IP + 1) % H2;
} else
last_i++;
buff_sinc = sinc_matrix[last_i];
chan_val = Complex(0.0, 0.0);
for (j=0; j<H2; j++)
chan_val += buff_f[(IP + j) % H2] * buff_sinc[j];
csi[t] = chan_val;
outp[t] = inp[t] * chan_val;
}
}
/***** Overloaded Function pass_through, providing AMPLITUDE CSI only *****/
void flat_rayleigh::pass_through(int length, Complex *inp,
Complex *outp, double *amp_csi)
{
register int k, j, t;
if (IndepFlag) {
chan_seed += 1; // reset the seed
// clean-up and then run idle, to ensure transients are dead
for (k=0; k<K; k++) {
st[0][k] = Complex(0.0, 0.0);
st[1][k] = Complex(0.0, 0.0);
st[2][k] = Complex(0.0, 0.0);
}
for (j=0; j<H2; j++)
buff_f[j] = Complex(0.0, 0.0);
// and now run idle for a little...
for (t=0; t<I*4*(4*K+H2); t++) {
if (last_i == I-1) {
last_i = 0;
/* produce new element, and insert it into buff_f[IP] */
st[0][0] = Complex( Gaussian(&chan_seed), Gaussian(&chan_seed) );
st[0][0] = st[0][0]*Ao - st[1][0]*c[0] - st[2][0]*d[0];
for (k=1; k<K; k++)
st[0][k] = st[0][k-1] + st[1][k-1]*a[k-1] + st[2][k-1]*b[k-1] -
st[1][k]*c[k] - st[2][k]*d[k];
buff_f[IP] = st[0][K-1] + st[1][K-1]*a[K-1] + st[2][K-1]*b[K-1];
for (k=0; k<K; k++) {
st[2][k] = st[1][k];
st[1][k] = st[0][k];
}
/* done, produced buff_f[IP], now move IP to next position */
IP = (IP + 1) % H2;
} else
last_i++;
// end here, don't produce chan_val, since idle run
}
} /* end of dry run for independence */
for (t=0; t<length; t++) {
if (last_i == I-1) {
last_i = 0;
/* now produce new element, and insert it into the buffer: buff_f[IP] */
st[0][0] = Complex( Gaussian(&chan_seed), Gaussian(&chan_seed) )*Ao -
st[1][0]*c[0] - st[2][0]*d[0];
for (k=1; k<K; k++)
st[0][k] = st[0][k-1] + st[1][k-1]*a[k-1] + st[2][k-1]*b[k-1] -
st[1][k]*c[k] - st[2][k]*d[k];
buff_f[IP] = st[0][K-1] + st[1][K-1]*a[K-1] + st[2][K-1]*b[K-1];
for (k=0; k<K; k++) {
st[2][k] = st[1][k];
st[1][k] = st[0][k];
}
/* done, produced buff_f[IP], now move IP to next position */
IP = (IP + 1) % H2;
} else
last_i++;
buff_sinc = sinc_matrix[last_i];
chan_val = Complex(0.0, 0.0);
for (j=0; j<H2; j++)
chan_val += buff_f[(IP + j) % H2] * buff_sinc[j];
amp_csi[t] = sqrt(abs(chan_val));
outp[t] = inp[t] * chan_val;
}
}
?? 快捷鍵說明
復制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -