解析matlab函数wrcoef的内部实现

启发：http://blog.sina.com.cn/s/blog_8ec096580101fsmn.html
一、首先，分析一个简单的例子，即对长度为2000的信号进行db6小波1层分解。
即：[C,L] = wavedec(s,1,'db6');A1 = wrcoef('a',C,L,'db6',1);
下面对函数wrcoef的实现过程进行详细解读：

第一步：将数组C中的CA1取出来，进行补零；

CA1 = C(1:L(1));
ca1 = zeros(1,2*L(1));
ca1(2:2:2*L(1)) = CA1(:);　
第二步：将CA1通过低通滤波器LO_R；
A11 = conv(LO_R,ca1);
第三步：得到1层分解后重构的概貌系数A1；
A1 = A11(12:2000+12-1);
注明：至于这里为什么要从第12个数据开始取，暂时不清楚。
从上面的分析来看，将分解最底层的概貌系数通过补零和滤波器的方式，直到得到和原始信号长度相同的数据。
二、分析3层分解的具体重构情况。Matlab语句：[C,L] = wavedec(s,3,'db6');
该句的作用是将信号s采用db6小波基进行3层分解，C里面中存放着小波变换后的各个分量的系数，L也是个数组，记录着C中每个元素的长度。
C = [CA3 CD3 CD2 CD1];L=[la3 ld3 ld2 ld1 N];(这里 L=[259 259 508 1005 2000]);
Matlab语句A3 = wrcoef('a',C,L,'db6',3);
该句的作用是利用C,L等数据信息，单支重构信号的概貌系数，得到原信号在该系数对应的尺度下的信号分量，其长度与原信号一致。
至于函数wrcoef的内部实现过程，可以看成是函数的wavedec的逆操作；下面对函数wrcoef的实现过程进行详细解读：

第一步：将数组C中的CA3取出来，长度为259，进行补零；
CA3= C(1:L(1));
ca3 = zeros(1,2*L1(1));
ca3(2:2:2*L1(1)) = CA3(:);
注明：
(1)wavedec在分解中采用偶数位置的欠采样，所以，这里是要在奇数位置补零；
(2)Wavedec函数中提到：假若滤波器的长度是2N，且信号长度是n，那么通过滤波器后的数据长度是(n+2N-1)，欠采样后得到的数据长度是floor((n-1)/2)+N
(3)补零前后的长度不是单纯的2倍关系。正如L中的数值所示，从开始的2000到1005，然后又从1005到508，依次类推。

第二步：长度259的数据通过补零+滤波器，得到长度为508的数据；
CA2 = conv(LO_R,CA33);
CA2 = CA2(12:519); % 519 = L(3)+12-1；
注明：wavedec在分解过程中用到的是LO_D，所以这里要使用对应的LO_R。

第三步：继续将得到的数据进行补零+滤波器(L(3)->L(4))；
ca2 = CA2(:);
CA2= zeros(1,2*L(3));
CA2(2:2:2* L(3)) = ca2(:);
CA1 = conv(LO_R,CA2);
CA1 = CA1(12:1016); % 1016=L(4)+12-1;

第四步：长度为L(4)的数据通过补零+滤波器，得到和信号长度相同的数据；
ca1 = CA1(:);
CA1 = zeros(1,2*L(4));
CA11(2:2:2*L(4)) = ca1(:);
A33 = conv(LO_R,CA1);
A33 = A33(12:2011); % 2011=L(5)+12-1;

第五步：比较推导得到的A33和直接公式得到的A3；
如果滤波器通过公式获取：[LO_D,HI_D,LO_R,HI_R] = wfilters('db6')，则
二者完全吻合；
如果滤波器如下定义：
LO_D = [ -0.0011 0.0048  0.0006  -0.0316 0.0275 0.0975 -0.1298 -0.2263 0.3153 0.7511 0.4946 0.1115];
LO_R = [ 0.1115  0.4946 0.7511  0.3153  -0.2263 -0.1298 0.0975  0.0275 -0.0316 0.0006 0.0048 -0.0011];
则数据A33和A3在数据的开头和结尾大约20个值误差较大。