mirror of
https://github.com/cookiengineer/audacity
synced 2025-08-16 08:34:10 +02:00
298 lines
6.9 KiB
C
298 lines
6.9 KiB
C
/* a few routines from a vector/matrix library */
|
|
#include "matlib.h"
|
|
|
|
void xpose(float *indata, long iRsiz, float *outdata, long oRsiz, long Nrows, long Ncols){
|
|
/* not in-place matrix transpose */
|
|
/* INPUTS */
|
|
/* *indata = input data array */
|
|
/* iRsiz = offset to between rows of input data array */
|
|
/* oRsiz = offset to between rows of output data array */
|
|
/* Nrows = number of rows in input data array */
|
|
/* Ncols = number of columns in input data array */
|
|
/* OUTPUTS */
|
|
/* *outdata = output data array */
|
|
|
|
float *irow; /* pointer to input row start */
|
|
float *ocol; /* pointer to output col start */
|
|
float *idata; /* pointer to input data */
|
|
float *odata; /* pointer to output data */
|
|
long RowCnt; /* row counter */
|
|
long ColCnt; /* col counter */
|
|
float T0; /* data storage */
|
|
float T1; /* data storage */
|
|
float T2; /* data storage */
|
|
float T3; /* data storage */
|
|
float T4; /* data storage */
|
|
float T5; /* data storage */
|
|
float T6; /* data storage */
|
|
float T7; /* data storage */
|
|
const long inRsizd1 = iRsiz;
|
|
const long inRsizd2 = 2*iRsiz;
|
|
const long inRsizd3 = inRsizd2+iRsiz;
|
|
const long inRsizd4 = 4*iRsiz;
|
|
const long inRsizd5 = inRsizd3+inRsizd2;
|
|
const long inRsizd6 = inRsizd4+inRsizd2;
|
|
const long inRsizd7 = inRsizd4+inRsizd3;
|
|
const long inRsizd8 = 8*iRsiz;
|
|
|
|
ocol = outdata;
|
|
irow = indata;
|
|
for (RowCnt=Nrows/8; RowCnt>0; RowCnt--){
|
|
idata = irow;
|
|
odata = ocol;
|
|
for (ColCnt=Ncols; ColCnt>0; ColCnt--){
|
|
T0 = *idata;
|
|
T1 = *(idata+inRsizd1);
|
|
T2 = *(idata+inRsizd2);
|
|
T3 = *(idata+inRsizd3);
|
|
T4 = *(idata+inRsizd4);
|
|
T5 = *(idata+inRsizd5);
|
|
T6 = *(idata+inRsizd6);
|
|
T7 = *(idata+inRsizd7);
|
|
*odata = T0;
|
|
*(odata+1) = T1;
|
|
*(odata+2) = T2;
|
|
*(odata+3) = T3;
|
|
*(odata+4) = T4;
|
|
*(odata+5) = T5;
|
|
*(odata+6) = T6;
|
|
*(odata+7) = T7;
|
|
idata++;
|
|
odata += oRsiz;
|
|
}
|
|
irow += inRsizd8;
|
|
ocol += 8;
|
|
}
|
|
if (Nrows%8 != 0){
|
|
for (ColCnt=Ncols; ColCnt>0; ColCnt--){
|
|
idata = irow++;
|
|
odata = ocol;
|
|
ocol += oRsiz;
|
|
for (RowCnt=Nrows%8; RowCnt>0; RowCnt--){
|
|
T0 = *idata;
|
|
*odata++ = T0;
|
|
idata += iRsiz;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
void cxpose(float *indata, long iRsiz, float *outdata, long oRsiz, long Nrows, long Ncols){
|
|
/* not in-place complex float matrix transpose */
|
|
/* INPUTS */
|
|
/* *indata = input data array */
|
|
/* iRsiz = offset to between rows of input data array */
|
|
/* oRsiz = offset to between rows of output data array */
|
|
/* Nrows = number of rows in input data array */
|
|
/* Ncols = number of columns in input data array */
|
|
/* OUTPUTS */
|
|
/* *outdata = output data array */
|
|
|
|
float *irow; /* pointer to input row start */
|
|
float *ocol; /* pointer to output col start */
|
|
float *idata; /* pointer to input data */
|
|
float *odata; /* pointer to output data */
|
|
long RowCnt; /* row counter */
|
|
long ColCnt; /* col counter */
|
|
float T0r; /* data storage */
|
|
float T0i; /* data storage */
|
|
float T1r; /* data storage */
|
|
float T1i; /* data storage */
|
|
float T2r; /* data storage */
|
|
float T2i; /* data storage */
|
|
float T3r; /* data storage */
|
|
float T3i; /* data storage */
|
|
const long inRsizd1 = 2*iRsiz;
|
|
const long inRsizd1i = 2*iRsiz + 1;
|
|
const long inRsizd2 = 4*iRsiz;
|
|
const long inRsizd2i = 4*iRsiz + 1;
|
|
const long inRsizd3 = inRsizd2+inRsizd1;
|
|
const long inRsizd3i = inRsizd2+inRsizd1 + 1;
|
|
const long inRsizd4 = 8*iRsiz;
|
|
|
|
ocol = outdata;
|
|
irow = indata;
|
|
for (RowCnt=Nrows/4; RowCnt>0; RowCnt--){
|
|
idata = irow;
|
|
odata = ocol;
|
|
for (ColCnt=Ncols; ColCnt>0; ColCnt--){
|
|
T0r = *idata;
|
|
T0i = *(idata +1);
|
|
T1r = *(idata+inRsizd1);
|
|
T1i = *(idata+inRsizd1i);
|
|
T2r = *(idata+inRsizd2);
|
|
T2i = *(idata+inRsizd2i);
|
|
T3r = *(idata+inRsizd3);
|
|
T3i = *(idata+inRsizd3i);
|
|
*odata = T0r;
|
|
*(odata+1) = T0i;
|
|
*(odata+2) = T1r;
|
|
*(odata+3) = T1i;
|
|
*(odata+4) = T2r;
|
|
*(odata+5) = T2i;
|
|
*(odata+6) = T3r;
|
|
*(odata+7) = T3i;
|
|
idata+=2;
|
|
odata += 2*oRsiz;
|
|
}
|
|
irow += inRsizd4;
|
|
ocol += 8;
|
|
}
|
|
if (Nrows%4 != 0){
|
|
for (ColCnt=Ncols; ColCnt>0; ColCnt--){
|
|
idata = irow;
|
|
odata = ocol;
|
|
for (RowCnt=Nrows%4; RowCnt>0; RowCnt--){
|
|
T0r = *idata;
|
|
T0i = *(idata+1);
|
|
*odata = T0r;
|
|
*(odata+1) = T0i;
|
|
odata+=2;
|
|
idata += 2*iRsiz;
|
|
}
|
|
irow+=2;
|
|
ocol += 2*oRsiz;
|
|
}
|
|
}
|
|
}
|
|
|
|
void cvprod(float *a, float *b, float *out, long N){
|
|
/* complex vector product, can be in-place */
|
|
/* product of complex vector *a times complex vector *b */
|
|
/* INPUTS */
|
|
/* N vector length */
|
|
/* *a complex vector length N complex numbers */
|
|
/* *b complex vector length N complex numbers */
|
|
/* OUTPUTS */
|
|
/* *out complex vector length N */
|
|
|
|
long OutCnt; /* output counter */
|
|
float A0R; /* A storage */
|
|
float A0I; /* A storage */
|
|
float A1R; /* A storage */
|
|
float A1I; /* A storage */
|
|
float A2R; /* A storage */
|
|
float A2I; /* A storage */
|
|
float A3R; /* A storage */
|
|
float A3I; /* A storage */
|
|
float B0R; /* B storage */
|
|
float B0I; /* B storage */
|
|
float B1R; /* B storage */
|
|
float B1I; /* B storage */
|
|
float B2R; /* B storage */
|
|
float B2I; /* B storage */
|
|
float B3R; /* B storage */
|
|
float B3I; /* B storage */
|
|
float T0R; /* TMP storage */
|
|
float T0I; /* TMP storage */
|
|
float T1R; /* TMP storage */
|
|
float T1I; /* TMP storage */
|
|
float T2R; /* TMP storage */
|
|
float T2I; /* TMP storage */
|
|
float T3R; /* TMP storage */
|
|
float T3I; /* TMP storage */
|
|
|
|
if (N>=4){
|
|
A0R = *a;
|
|
B0R = *b;
|
|
A0I = *(a +1);
|
|
B0I = *(b +1);
|
|
A1R = *(a +2);
|
|
B1R = *(b +2);
|
|
A1I = *(a +3);
|
|
B1I = *(b +3);
|
|
A2R = *(a +4);
|
|
B2R = *(b +4);
|
|
A2I = *(a +5);
|
|
B2I = *(b +5);
|
|
A3R = *(a +6);
|
|
B3R = *(b +6);
|
|
A3I = *(a +7);
|
|
B3I = *(b +7);
|
|
T0R = A0R * B0R;
|
|
T0I = (A0R * B0I);
|
|
T1R = A1R * B1R;
|
|
T1I = (A1R * B1I);
|
|
T2R = A2R * B2R;
|
|
T2I = (A2R * B2I);
|
|
T3R = A3R * B3R;
|
|
T3I = (A3R * B3I);
|
|
T0R -= (A0I * B0I);
|
|
T0I = A0I * B0R + T0I;
|
|
T1R -= (A1I * B1I);
|
|
T1I = A1I * B1R + T1I;
|
|
T2R -= (A2I * B2I);
|
|
T2I = A2I * B2R + T2I;
|
|
T3R -= (A3I * B3I);
|
|
T3I = A3I * B3R + T3I;
|
|
for (OutCnt=N/4-1; OutCnt > 0; OutCnt--){
|
|
a += 8;
|
|
b += 8;
|
|
A0R = *a;
|
|
B0R = *b;
|
|
A0I = *(a +1);
|
|
B0I = *(b +1);
|
|
A1R = *(a +2);
|
|
B1R = *(b +2);
|
|
A1I = *(a +3);
|
|
B1I = *(b +3);
|
|
A2R = *(a +4);
|
|
B2R = *(b +4);
|
|
A2I = *(a +5);
|
|
B2I = *(b +5);
|
|
A3R = *(a +6);
|
|
B3R = *(b +6);
|
|
A3I = *(a +7);
|
|
B3I = *(b +7);
|
|
*out = T0R;
|
|
*(out +1) = T0I;
|
|
*(out +2) = T1R;
|
|
*(out +3) = T1I;
|
|
*(out +4) = T2R;
|
|
*(out +5) = T2I;
|
|
*(out +6) = T3R;
|
|
*(out +7) = T3I;
|
|
T0R = A0R * B0R;
|
|
T0I = (A0R * B0I);
|
|
T1R = A1R * B1R;
|
|
T1I = (A1R * B1I);
|
|
T2R = A2R * B2R;
|
|
T2I = (A2R * B2I);
|
|
T3R = A3R * B3R;
|
|
T3I = (A3R * B3I);
|
|
T0R -= (A0I * B0I);
|
|
T0I = A0I * B0R + T0I;
|
|
T1R -= (A1I * B1I);
|
|
T1I = A1I * B1R + T1I;
|
|
T2R -= (A2I * B2I);
|
|
T2I = A2I * B2R + T2I;
|
|
T3R -= (A3I * B3I);
|
|
T3I = A3I * B3R + T3I;
|
|
out += 8;
|
|
}
|
|
a += 8;
|
|
b += 8;
|
|
*out = T0R;
|
|
*(out +1) = T0I;
|
|
*(out +2) = T1R;
|
|
*(out +3) = T1I;
|
|
*(out +4) = T2R;
|
|
*(out +5) = T2I;
|
|
*(out +6) = T3R;
|
|
*(out +7) = T3I;
|
|
out += 8;
|
|
}
|
|
for (OutCnt=N%4; OutCnt > 0; OutCnt--){
|
|
A0R = *a++;
|
|
B0R = *b++;
|
|
A0I = *a++;
|
|
B0I = *b++;
|
|
T0R = A0R * B0R;
|
|
T0I = (A0R * B0I);
|
|
T0R -= (A0I * B0I);
|
|
T0I = A0I * B0R + T0I;
|
|
*out++ = T0R;
|
|
*out++ = T0I;
|
|
}
|
|
}
|