?? dwt.c
字號:
/*--------------------------------------------------------------------------------------------*/
#include "dwt.h"
/*--------------------------------------------------------------------------------------------*/
static REAL *LineA, *LineB, *LineC;
/* The synthesis delay is at most half of analysis delay - upsampling */
#define MaxDelay 11
/* Note that the W53 and W137 are not normalized, so their uses in bitplane coder
* need to be compensated with weighting. In fact, it is better to use the
* normalized coefficients directly.
*/
#define D97 0
#define V1018 1
#define ADEL 2
#define W53 3
#define W137 4
char *SubbandFilterName[5] = {"antonini", "villa1810", "ADEL", "W5x3", "W13x7"};
/*--------------------------------------------------------------------------------------------*/
/* Daubechies 9/7 wavelets */
#if 1
static REAL AL97[5] = {
(REAL)8.526986790094022e-01,
(REAL)3.774028556126537e-01,
(REAL)-1.106244044184226e-01,
(REAL)-2.384946501937986e-02,
(REAL)3.782845550699535e-02
};
static REAL AH97[4] = {
(REAL)7.884856164056651e-01,
(REAL)-4.180922732222124e-01,
(REAL)-4.068941760955867e-02,
(REAL)6.453888262893856e-02
};
static REAL SL97[4] = {
(REAL)7.884856164056651e-01,
(REAL)4.180922732222124e-01,
(REAL)-4.068941760955867e-02,
};
static REAL SH97[5] = {
(REAL)8.526986790094022e-01,
(REAL)-3.774028556126537e-01,
(REAL)-1.106244044184226e-01,
(REAL)2.384946501937986e-02,
(REAL)3.782845550699535e-02
};
#else
/*-------------------------------------------------------*/
static REAL AL97[5] = {
(REAL)(8.526986790094022e-01/M_SQRT2),
(REAL)(3.774028556126537e-01/M_SQRT2),
(REAL)(-1.106244044184226e-01/M_SQRT2),
(REAL)(-2.384946501937986e-02/M_SQRT2),
(REAL)(3.782845550699535e-02/M_SQRT2)
};
static REAL AH97[4] = {
(REAL)(7.884856164056651e-01/M_SQRT2),
(REAL)(-4.180922732222124e-01/M_SQRT2),
(REAL)(-4.068941760955867e-02/M_SQRT2),
(REAL)(6.453888262893856e-02/M_SQRT2)
};
static REAL SL97[4] = {
(REAL)(7.884856164056651e-01/M_SQRT2),
(REAL)(4.180922732222124e-01/M_SQRT2),
(REAL)(-4.068941760955867e-02/M_SQRT2),
(REAL)(-6.453888262893856e-02/M_SQRT2)
};
static REAL SH97[5] = {
(REAL)(8.526986790094022e-01/M_SQRT2),
(REAL)(-3.774028556126537e-01/M_SQRT2),
(REAL)(-1.106244044184226e-01/M_SQRT2),
(REAL)(2.384946501937986e-02/M_SQRT2),
(REAL)(3.782845550699535e-02/M_SQRT2)
};
#endif
/*--------------------------------------------------------------------------------------------*/
/* Le-Gall 5/3 */
/* The rational coeff : -1/8 1/4 3/4 1/4 -1/8 */
/* normalized by *sqrt(2) is as shown below */
//static REAL AL53[3] = { 1.060660171779821e+00, 3.535533905932738e-01, -1.767766952966369e-01};
static REAL AL53[3] = { 0.75, 0.25, -0.125};
//static REAL AH53[2] = { 7.071067811865476e-01, -3.535533905932738e-01};
static REAL AH53[2] = { 1, -0.5};
/* rational coeff: 1/2 1 1/2 */
/* normalized by 1/sqrt(2) as shown below */
//static REAL SL53[2] = { 7.071067811865476e-01, 3.535533905932738e-01};
static REAL SL53[2] = { 1, 0.5};
//static REAL SH53[3] = { 1.060660171779821e+00, -3.535533905932738e-01, -1.767766952966369e-01};
static REAL SH53[3] = { 0.75, -0.25, -0.125};
/*--------------------------------------------------------------------------------------------*/
static REAL AL137[7] = {0.6796875, 0.28125, -0.123046875, -0.03125, 0.03515625, 0, -1.953125e-3};
static REAL AH137[4] = {1, -0.5625, 0, 0.0625};
static REAL SL137[4] = {1, 0.5625, 0, -0.0625};
static REAL SH137[7] = {0.6796875, -0.28125, -0.123046875, 0.03125, 0.03515625, 0, -1.953125e-3};
/*--------------------------------------------------------------------------------------------*/
/* Villasenor 10/18 wavelets */
static REAL AL1018[5] = {
(REAL)7.589077294537619e-01,
(REAL)7.679048884691436e-02,
(REAL)-1.575264469076351e-01,
(REAL)8.244478227504624e-05,
(REAL)2.885256501123136e-02
};
static REAL AH1018 [9] = {
(REAL)-6.233596410344158e-01,
(REAL)1.633685405569888e-01,
(REAL)8.566118833165885e-02,
(REAL)-1.376513483818652e-02,
(REAL)-3.083373438534267e-02,
(REAL)-2.528037293949898e-03,
(REAL)9.452462998353147e-03,
(REAL)-2.727196296995984e-06,
(REAL)-9.544158682436510e-04
};
static REAL SL1018 [9] = {
(REAL)6.233596410344158e-01,
(REAL)1.633685405569888e-01,
(REAL)-8.566118833165885e-02,
(REAL)-1.376513483818652e-02,
(REAL)3.083373438534267e-02,
(REAL)-2.528037293949898e-03,
(REAL)-9.452462998353147e-03,
(REAL)-2.727196296995984e-06,
(REAL)9.544158682436510e-04
};
static REAL SH1018[5] = {
(REAL)7.589077294537619e-01,
(REAL) -7.679048884691436e-02,
(REAL)-1.575264469076351e-01,
(REAL)-8.244478227504624e-05,
(REAL)2.885256501123136e-02
};
/*--------------------------------------------------------------------------------------------*/
#define _EPIC_QMF_ 1
#if _EPIC_QMF_
static REAL AdelAL[5] = {
(REAL) 0.7973282005652047,
(REAL) 0.4146915396805326,
(REAL)-0.07338062349083309,
(REAL)-0.06093975981002458,
(REAL) 0.028071524524270521
};
static REAL AdelAH[5] = {
(REAL) 0.7973282005652047,
(REAL)-0.4146915396805326,
(REAL)-0.07338062349083309,
(REAL) 0.06093975981002458,
(REAL) 0.028071524524270521
};
#else
static REAL AdelAL[5] = {
(REAL) 0.798436693,
(REAL) 0.413955551,
(REAL)-0.073878516,
(REAL)-0.060401061,
(REAL) 0.02821356
};
static REAL AdelAH[5] = {
(REAL) 0.798436693,
(REAL)-0.413955551,
(REAL)-0.073878516,
(REAL) 0.060401061,
(REAL) 0.02821356
};
#endif
static REAL *AdelSL=AdelAL;
static REAL *AdelSH=AdelAH;
/*--------------------------------------------------------------------------------------------*/
/* From A.Said & W. Pearlman SPIHT code */
int SubbandMaxScales(int n)
{
int l1, l2;
for (l1 = 0; !(n & 1); l1++) n >>= 1;
for (l2 = l1 - 4; n; l2++) n >>= 1;
return (l1 < l2 ? l1 : l2);
}
/*--------------------------------------------------------------------------------------------*/
void InitializeSubbandBuffers(int rows, int cols)
{
int n;
/* allocate the larger of row or col */
if (rows>cols){
n = rows;
}
else{
n = cols;
}
/* temporary buffer for input and output of transform data */
/* we are allocating more than enough memory */
/*
* For the 2-D transform, we have the following settings:
*
* Analysis:
* input - LineA zero offset data.
* output - LineB start zero offset; first half low-pass, second half high-pass.
*
* Synthesis:
* col: input - LineA zero offset data.
* output - LineA zero offset.
* row: input - current row with proper size.
* output - current row with proper size.
*/
/* For 1-D transform, the input to analysis and synthesis filters can be
* any 1-D array. The output of analysis filterings are in LineB.
* The output of the synthesis is the input array itself.
*/
if ((LineA = (REAL *)malloc(sizeof(REAL)*(n+2*MaxDelay)))==NULL){
Error("Fail to allocate memory.\n");
}
if ((LineB = (REAL *)malloc(sizeof(REAL)*(n+2*MaxDelay)))==NULL){
Error("Fail to allocate memory.\n");
}
if ((LineC = (REAL *)malloc(sizeof(REAL)*(n+2*MaxDelay)))==NULL){
Error("Fail to allocate memory.\n");
}
}
/*--------------------------------------------------------------------------------------------*/
void FreeSubbandBuffers(void)
{
free(LineA);
free(LineB);
free(LineC);
}
/*--------------------------------------------------------------------------------------------*/
#define AnalDelay97 4
#define SyntDelay97 2
/* 1-D analysis filter bank - odd length LP FIR filter */
void SubbandAnal97(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 + AnalDelay97 + 1;
LinePtrL = InPtrL - 2;
LinePtrR = in + n + AnalDelay97;
InPtrR = LinePtrR - 2;
for (i=0; i<AnalDelay97; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
halfL = (n+1)>>1;
/* low-pass */
for (i=0, j=AnalDelay97; i<halfL; i++, j+=2){
LineB[i] = AL97[0]*in[j]
+ AL97[1]*(in[j-1]+in[j+1])
+ AL97[2]*(in[j-2]+in[j+2])
+ AL97[3]*(in[j-3]+in[j+3])
+ AL97[4]*(in[j-4]+in[j+4]);
}
/* high-pass - align the filter to odd index samples */
for (i=halfL, j=AnalDelay97+1; i<n; i++, j+=2){
LineB[i] = AH97[0]*in[j]
+ AH97[1]*(in[j-1]+in[j+1])
+ AH97[2]*(in[j-2]+in[j+2])
+ AH97[3]*(in[j-3]+in[j+3]);
}
}
/*--------------------------------------------------------------------------------------------*/
/* n is the total signal length.
* in contains the data (zero offset).
*/
void SubbandSynt97(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+SyntDelay97;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelay97+1;
LinePtrL = InPtrL-2;
InPtrR = LineB+SyntDelay97+halfL-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelay97; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (H, H) */
/* copy the subband */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelay97;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelay97;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelay97+halfH-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelay97; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
}
else{
/* Low-pass: (W, H) */
/* copy the subband */
InPtrL = in;
LinePtrL = LineB+SyntDelay97;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelay97+1;
LinePtrL = InPtrL-2;
InPtrR = LineB+SyntDelay97+halfL-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelay97; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (H, W) */
/* copy the subband */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelay97;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelay97;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelay97+halfH-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelay97; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
}
for (i=0, j=SyntDelay97; j<halfH+SyntDelay97; i+=2, j++){
in[i] = SL97[0] * LineB[j]
+ SL97[2] * (LineB[j+1] + LineB[j-1])
+ SH97[1] * (LineC[j] + LineC[j-1])
+ SH97[3] * (LineC[j+1] + LineC[j-2]);
in[i+1] = SL97[1] * (LineB[j] + LineB[j+1])
+ SL97[3] * (LineB[j-1] + LineB[j+2])
+ SH97[0] * LineC[j]
+ SH97[2] * (LineC[j+1] + LineC[j-1])
+ SH97[4] * (LineC[j+2] + LineC[j-2]);
}
if (Odd){
i = halfH<<1;
j = halfH+SyntDelay97;
in[i] = SL97[0] * LineB[j]
+ SL97[2] * (LineB[j+1] + LineB[j-1])
+ SH97[1] * (LineC[j] + LineC[j-1])
+ SH97[3] * (LineC[j+1] + LineC[j-2]);
}
}
/*--------------------------------------------------------------------------------------------*/
#undef AnalDelay97
#undef SyntDelay97
/*--------------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------------*/
#define AnalDelay137 6
#define SyntDelay137 3
/*--------------------------------------------------------------------------------------------*/
/* 1-D analysis filter bank - odd length LP FIR filter */
void SubbandAnal137(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 + AnalDelay137 + 1;
LinePtrL = InPtrL - 2;
LinePtrR = in + n + AnalDelay137;
InPtrR = LinePtrR - 2;
for (i=0; i<AnalDelay137; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
halfL = (n+1)>>1;
/* low-pass */
for (i=0, j=AnalDelay137; i<halfL; i++, j+=2){
LineB[i] = AL137[0]*in[j]
+ AL137[1]*(in[j-1]+in[j+1])
+ AL137[2]*(in[j-2]+in[j+2])
+ AL137[3]*(in[j-3]+in[j+3])
+ AL137[4]*(in[j-4]+in[j+4])
+ AL137[5]*(in[j-5]+in[j+5])
+ AL137[6]*(in[j-6]+in[j+6]);
}
/* high-pass - align the filter to odd index samples */
for (i=halfL, j=AnalDelay137+1; i<n; i++, j+=2){
LineB[i] = AH137[0]*in[j]
+ AH137[1]*(in[j-1]+in[j+1])
+ AH137[2]*(in[j-2]+in[j+2])
+ AH137[3]*(in[j-3]+in[j+3]);
}
}
/*--------------------------------------------------------------------------------------------*/
/* n is the total signal length.
* in contains the data (zero offset).
*/
void SubbandSynt137(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+SyntDelay137;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelay137+1;
LinePtrL = InPtrL-2;
InPtrR = LineB+SyntDelay137+halfL-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelay137; i++){
*LinePtrL-- = *InPtrL++;
}
/* High-pass: (H, H) */
/* copy the subband */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelay137;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelay137;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelay137+halfH-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelay137; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
}
else{
/* Low-pass: (W, H) */
/* copy the subband */
InPtrL = in;
LinePtrL = LineB+SyntDelay137;
for (i=0; i<halfL; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineB+SyntDelay137+1;
LinePtrL = InPtrL-2;
InPtrR = LineB+SyntDelay137+halfL-1;
LinePtrR = InPtrR+1;
for (i=0; i<SyntDelay137; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
}
/* High-pass: (H, W) */
/* copy the subband */
InPtrL = in+halfL;
LinePtrL = LineC+SyntDelay137;
for (i=0; i<halfH; i++){
*LinePtrL++ = *InPtrL++;
}
InPtrL = LineC+SyntDelay137;
LinePtrL = InPtrL-1;
InPtrR = LineC+SyntDelay137+halfH-2;
LinePtrR = InPtrR+2;
for (i=0; i<SyntDelay137; i++){
*LinePtrL-- = *InPtrL++;
*LinePtrR++ = *InPtrR--;
?? 快捷鍵說明
復(fù)制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -