📄 dwt.c
字号:
}
}
for (i=0, j=SyntDelay137; j<halfH+SyntDelay137; i+=2, j++){
in[i] = SL137[0] * LineB[j]
+ SL137[2] * (LineB[j+1] + LineB[j-1])
+ SH137[1] * (LineC[j] + LineC[j-1])
+ SH137[3] * (LineC[j+1] + LineC[j-2]);
in[i+1] = SL137[1] * (LineB[j] + LineB[j+1])
+ SL137[3] * (LineB[j-1] + LineB[j+2])
+ SH137[0] * LineC[j]
+ SH137[2] * (LineC[j+1] + LineC[j-1])
+ SH137[4] * (LineC[j+2] + LineC[j-2])
+ SH137[6] * (LineC[j+3] + LineC[j-3]);
}
if (Odd){
i = halfH<<1;
j = halfH+SyntDelay137;
in[i] = SL137[0] * LineB[j]
+ SL137[2] * (LineB[j+1] + LineB[j-1])
+ SH137[1] * (LineC[j] + LineC[j-1])
+ SH137[3] * (LineC[j+1] + LineC[j-2]);
}
}
/*--------------------------------------------------------------------------------------------*/
#undef AnalDelay137
#undef SyntDelay137
/*--------------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------------*/
#define AnalDelay53 2
#define SyntDelay53 1
/*--------------------------------------------------------------------------------------------*/
void SubbandAnal53(REAL *in, int n)
{
int i, j, halfL;
REAL *LinePtr, *InPtr, *LinePtrL, *LinePtrR, *InPtrL, *InPtrR;
/* n is the data size (not including extension) */
/* data in n start at AnalDelay */
/* (W, W) */
InPtrL = in + AnalDelay53 + 1;
LinePtrL = InPtrL - 2;
LinePtrR = in + n + AnalDelay53;
InPtrR = LinePtrR - 2;
for (i=0; i<AnalDelay53; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
halfL = (n+1)>>1;
/* low-pass */
for (i=0, j=AnalDelay53; i<halfL; i++, j+=2){
LineB[i] = AL53[0]*in[j]
+ AL53[1]*(in[j-1]+in[j+1])
+ AL53[2]*(in[j-2]+in[j+2]);
}
/* high-pass - align the filter to odd index samples */
for (i=halfL, j=AnalDelay53+1; i<n; i++, j+=2){
LineB[i] = AH53[0]*in[j]
+ AH53[1]*(in[j-1]+in[j+1]);
}
}
/*--------------------------------------------------------------------------------------------*/
void SubbandSynt53(REAL *in, int n)
{
int i, j, halfL, halfH;
int Odd;
REAL *LinePtr, *InPtr, *LinePtrL, *LinePtrR, *InPtrL, *InPtrR;
halfL = (n+1)>>1;
halfH = (n)>>1;
Odd = (n & 0x1);
if (Odd){
/* Low-pass: (W, W) */
/* copy the subband */
InPtrL = in;
LinePtrL = LineB+SyntDelay53;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelay53+1;
LinePtrL = InPtrL-2;
InPtrR = LineB+SyntDelay53+halfL-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelay53; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (H, H) */
/* copy the subband */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelay53;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelay53;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelay53+halfH-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelay53; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
}
else{
/* Low-pass: (W, H) */
/* copy the subband */
InPtrL = in;
LinePtrL = LineB+SyntDelay53;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelay53+1;
LinePtrL = InPtrL-2;
InPtrR = LineB+SyntDelay53+halfL-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelay53; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (H, W) */
/* copy the subband */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelay53;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelay53;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelay53+halfH-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelay53; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
}
for (i=0, j=SyntDelay53; j<halfH+SyntDelay53; i+=2, j++){
in[i] = SL53[0] * LineB[j]
+ SH53[1] * (LineC[j] + LineC[j-1]);
in[i+1] = SL53[1] * (LineB[j] + LineB[j+1])
+ SH53[0] * LineC[j]
+ SH53[2] * (LineC[j+1] + LineC[j-1]);
}
if (Odd){
i = halfH<<1;
j = halfH+SyntDelay53;
in[i] = SL53[0] * LineB[j]
+ SH53[1] * (LineC[j] + LineC[j-1]);
}
}
/*--------------------------------------------------------------------------------------------*/
#undef AnalDelay53
#undef SyntDelay53
/*--------------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------------*/
#define AnalDelay1018 8
#define SyntDelay1018 4
/*--------------------------------------------------------------------------------------------*/
/* 1-D analysis filter bank - even length LP FIR filter */
void SubbandAnal1018(REAL *in, int n)
{
int i, j, halfL;
REAL *LinePtr, *InPtr, *LinePtrL, *LinePtrR, *InPtrL, *InPtrR;
/* n is the data size (not including extension) */
/* data in n start at AnalDelay */
/* (H, H) */
InPtrL = in + AnalDelay1018;
LinePtrL = InPtrL-1;
LinePtrR = in + n + AnalDelay1018;
InPtrR = LinePtrR-1;
for (i=0; i<AnalDelay1018; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
halfL = (n+1)>>1;
/* low-pass */
for (i=0, j=AnalDelay1018; i<halfL; i++, j+=2){
LineB[i] = AL1018[0]*(in[j] + in[j+1])
+ AL1018[1]*(in[j-1] + in[j+2])
+ AL1018[2]*(in[j-2] + in[j+3])
+ AL1018[3]*(in[j-3] + in[j+4])
+ AL1018[4]*(in[j-4] + in[j+5]);
}
/* high-pass */
for (i=halfL, j=AnalDelay1018; i<n; i++, j+=2){
LineB[i] = AH1018[0]*(in[j] - in[j+1])
+ AH1018[1]*(in[j-1] - in[j+2])
+ AH1018[2]*(in[j-2] - in[j+3])
+ AH1018[3]*(in[j-3] - in[j+4])
+ AH1018[4]*(in[j-4] - in[j+5])
+ AH1018[5]*(in[j-5] - in[j+6])
+ AH1018[6]*(in[j-6] - in[j+7])
+ AH1018[7]*(in[j-7] - in[j+8])
+ AH1018[8]*(in[j-8] - in[j+9]);
}
}
/*--------------------------------------------------------------------------------------------*/
/* n is the total signal length.
* in contains the data (zero offset).
*/
void SubbandSynt1018(REAL *in, int n)
{
int i, j, halfL, halfH;
int Odd;
REAL *LinePtr, *InPtr, *LinePtrL, *LinePtrR, *InPtrL, *InPtrR;
halfL = (n+1)>>1;
halfH = (n)>>1;
Odd = (n & 0x1);
if (Odd){
/* Low-pass: (H, W) */
/* copy the subband */
InPtrL = in;
LinePtrL = LineB+SyntDelay1018;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelay1018;
LinePtrL = InPtrL-1;
InPtrR = LineB+SyntDelay1018+halfL-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelay1018; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (-H, +[0]-W) */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelay1018;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
/* insert 0 */
*LinePtrL = 0;
InPtrL = LineC+SyntDelay1018;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelay1018+halfH-1;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelay1018; i++){
*LinePtrL-- = -(*InPtrL++);
*LinePtrR++ = -(*InPtrR--);
}
}
else{
/* Low-pass: (H, H) */
/* copy the subband */
InPtrL = in;
LinePtrL = LineB+SyntDelay1018;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelay1018;
LinePtrL = InPtrL-1;
InPtrR = LineB+SyntDelay1018+halfL-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelay1018; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (-H, -H) */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelay1018;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelay1018;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelay1018+halfH-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelay1018; i++){
*LinePtrL-- = -(*InPtrL++);
*LinePtrR++ = -(*InPtrR--);
}
}
/* there is no way to exploit the symmetry */
for (i=0, j=SyntDelay1018; j<halfH+SyntDelay1018; i+=2, j++){
in[i] = SL1018[0] * LineB[j]
+ SL1018[1] * LineB[j-1]
+ SL1018[3] * LineB[j-2]
+ SL1018[5] * LineB[j-3]
+ SL1018[7] * LineB[j-4]
+ SL1018[2] * LineB[j+1]
+ SL1018[4] * LineB[j+2]
+ SL1018[6] * LineB[j+3]
+ SL1018[8] * LineB[j+4]
- SH1018[0] * LineC[j]
+ SH1018[1] * LineC[j-1]
+ SH1018[3] * LineC[j-2]
- SH1018[2] * LineC[j+1]
- SH1018[4] * LineC[j+2];
in[i+1] = SL1018[0] * LineB[j]
+ SL1018[2] * LineB[j-1]
+ SL1018[4] * LineB[j-2]
+ SL1018[6] * LineB[j-3]
+ SL1018[8] * LineB[j-4]
+ SL1018[1] * LineB[j+1]
+ SL1018[3] * LineB[j+2]
+ SL1018[5] * LineB[j+3]
+ SL1018[7] * LineB[j+4]
+ SH1018[0] * LineC[j]
+ SH1018[2] * LineC[j-1]
+ SH1018[4] * LineC[j-2]
- SH1018[1] * LineC[j+1]
- SH1018[3] * LineC[j+2];
}
if (Odd){
i = halfH<<1;
j = halfH+SyntDelay1018;
in[i] = SL1018[0] * LineB[j]
+ SL1018[1] * LineB[j-1]
+ SL1018[3] * LineB[j-2]
+ SL1018[5] * LineB[j-3]
+ SL1018[7] * LineB[j-4]
+ SL1018[2] * LineB[j+1]
+ SL1018[4] * LineB[j+2]
+ SL1018[6] * LineB[j+3]
+ SL1018[8] * LineB[j+4]
- SH1018[0] * LineC[j]
+ SH1018[1] * LineC[j-1]
+ SH1018[3] * LineC[j-2]
- SH1018[2] * LineC[j+1]
- SH1018[4] * LineC[j+2];
}
}
/*--------------------------------------------------------------------------------------------*/
#undef AnalDelay1018
#undef SyntDelay1018
/*--------------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------------*/
#define AnalDelayAdel 4
#define SyntDelayAdel 4
/*--------------------------------------------------------------------------------------------*/
/* 1-D analysis filter bank - odd length QMF FIR filter */
void SubbandAnalAdel(REAL *in, int n)
{
int i, j, halfL;
REAL *LinePtr, *InPtr, *LinePtrL, *LinePtrR, *InPtrL, *InPtrR;
/* n is the data size (not including extension) */
/* data in n start at AnalDelay */
/* (W, W) */
InPtrL = in + AnalDelayAdel + 1;
LinePtrL = InPtrL - 2;
LinePtrR = in + n + AnalDelayAdel;
InPtrR = LinePtrR - 2;
for (i=0; i<AnalDelayAdel; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
halfL = (n+1)>>1;
/* low-pass */
for (i=0, j=AnalDelayAdel; i<halfL; i++, j+=2){
LineB[i] = AdelAL[0]*in[j]
+ AdelAL[1]*(in[j-1]+in[j+1])
+ AdelAL[2]*(in[j-2]+in[j+2])
+ AdelAL[3]*(in[j-3]+in[j+3])
+ AdelAL[4]*(in[j-4]+in[j+4]);
}
/* high-pass - align the filter to odd index samples */
for (i=halfL, j=AnalDelayAdel+1; i<n; i++, j+=2){
LineB[i] = AdelAH[0]*in[j]
+ AdelAH[1]*(in[j-1]+in[j+1])
+ AdelAH[2]*(in[j-2]+in[j+2])
+ AdelAH[3]*(in[j-3]+in[j+3])
+ AdelAH[4]*(in[j-4]+in[j+4]);
}
}
/*--------------------------------------------------------------------------------------------*/
/* n is the total signal length.
* in contains the data (zero offset).
*/
void SubbandSyntAdel(REAL *in, int n)
{
int i, j, halfL, halfH;
int Odd;
REAL *LinePtr, *InPtr, *LinePtrL, *LinePtrR, *InPtrL, *InPtrR;
halfL = (n+1)>>1;
halfH = (n)>>1;
Odd = (n & 0x1);
if (Odd){
/* Low-pass: (W, W) */
/* copy the subband */
InPtrL = in;
LinePtrL = LineB+SyntDelayAdel;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelayAdel+1;
LinePtrL = InPtrL-2;
InPtrR = LineB+SyntDelayAdel+halfL-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelayAdel; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (H, H) */
/* copy the subband */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelayAdel;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelayAdel;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelayAdel+halfH-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelayAdel; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
}
else{
/* Low-pass: (W, H) */
/* copy the subband */
InPtrL = in;
LinePtrL = LineB+SyntDelayAdel;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelayAdel+1;
LinePtrL = InPtrL-2;
InPtrR = LineB+SyntDelayAdel+halfL-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelayAdel; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (H, W) */
/* copy the subband */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelayAdel;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelayAdel;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelayAdel+halfH-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelayAdel; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
}
for (i=0, j=SyntDelayAdel; j<halfH+SyntDelayAdel; i+=2, j++){
in[i] = AdelSL[0] * LineB[j]
+ AdelSL[2] * (LineB[j+1] + LineB[j-1])
+ AdelSL[4] * (LineB[j+2] + LineB[j-2])
+ AdelSH[1] * (LineC[j] + LineC[j-1])
+ AdelSH[3] * (LineC[j+1] + LineC[j-2]);
in[i+1] = AdelSL[1] * (LineB[j] + LineB[j+1])
+ AdelSL[3] * (LineB[j-1] + LineB[j+2])
+ AdelSH[0] * LineC[j]
+ AdelSH[2] * (LineC[j+1] + LineC[j-1])
+ AdelSH[4] * (LineC[j+2] + LineC[j-2]);
}
if (Odd){
i = halfH<<1;
j = halfH+SyntDelayAdel;
in[i] = AdelSL[0] * LineB[j]
+ AdelSL[2] * (LineB[j+1] + LineB[j-1])
+ AdelSL[4] * (LineB[j+2] + LineB[j-2])
+ AdelSH[1] * (LineC[j] + LineC[j-1])
+ AdelSH[3] * (LineC[j+1] + LineC[j-2]);
}
}
/*--------------------------------------------------------------------------------------------*/
#undef AnalAdel
#undef SyntAdel
/*--------------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------------*/
DWTransform *DWTransformAlloc(int XSize, int YSize, int nLevels)
{
int i;
DWTransform *t;
int *LowXSize, *LowYSize, *HighXSize, *HighYSize;
if ((t=(DWTransform *)malloc(sizeof(DWTransform)))==NULL){
return NULL;
}
t->XSize = XSize;
t->YSize = YSize;
t->nLevels = nLevels;
t->nSubbands = 3*t->nLevels+1;
t->Coeff=(REAL *)calloc(t->XSize*t->YSize, sizeof(REAL));
t->SubbandSize = (int *)calloc(t->nSubbands, sizeof(int));
t->SubbandXSize = (int *)calloc(t->nSubbands, sizeof(int));
t->SubbandYSize = (int *)calloc(t->nSubbands, sizeof(int));
t->SubbandScale = (int *)calloc(t->nSubbands, sizeof(int));
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -