實驗3-率失真函數(shù)計算的程序設(shè)計

實驗3-率失真函數(shù)計算的程序設(shè)計實驗3-率失真函數(shù)計算的程序設(shè)計實驗3-率失真函數(shù)計算的程序設(shè)計xxx公司實驗3-率失真函數(shù)計算的程序設(shè)計文件編號：文件日期：修訂次數(shù)：第1.0次更改批準審核制定方案設(shè)計，管理制度實驗3率失真函數(shù)計算的程序設(shè)計一、實驗問題假定一個離散無記憶信源（DMS）的信源符號集為Au={1,2,...,r},其概率分布為p(u)；信宿符號集為Av={1,2,....,s}。而失真?zhèn)榷染仃嚍橐粋€rs維矩陣D=[dij]。利用Matlab畫出率失真函數(shù)R()的曲線圖。二、實驗環(huán)境計算機、Windows2000或以上、Matlab或以上三、實驗目的了解率失真函數(shù)性質(zhì)、意義。掌握簡單的率失真函數(shù)計算方法；掌握使用Matlab實現(xiàn)一般率失真函數(shù)的計算方法；掌握Matlab求解非線性方程組的方法。四、實驗內(nèi)容從理論上計算r=s=2。p(u=1)=p,p(u=2)=1-p;d=[0,1;1,0]的率失真函數(shù)R()。對一般性的DMS信源，計算率失真函數(shù)R()的理論公式進行推導。找出比較合適的方程求解方法。使用編制Matlab編制程序求解一般的率失真函數(shù)R()。給定r=s=2。p(u=1)=,p=(u=2)=;d=[0,1;1,0]，測試程序，即比較程序運行結(jié)果與理論計算結(jié)果,改變參數(shù)，畫出函數(shù)圖。顯示在計算精度為以及運行計算的配置(CPU型號、CPU的頻率、內(nèi)存的)的條件下，系統(tǒng)循環(huán)次數(shù)、累計計算時間、平均每次循環(huán)所用時間等。五、實驗要求提前預習實驗，認真閱讀實驗原理。認真高效的完成實驗，實驗過程中服從實驗室管理人員以及實驗指導老師的管理。將實驗報告寫成論文的形式。要求有：問題的提出：包括R()的物理意義、用途（可以舉出具體的用途）、計算的困難性等。解決問題的原理方法：包括所有的公式推導的細節(jié)。解決問題的具體方法：包括程序框圖及Matlab源程序。實驗結(jié)果：利用你的程序給出不同參數(shù)得到的實驗結(jié)果，包括實驗曲線圖、程序循環(huán)次數(shù)、累計計算時間、平均每次循環(huán)所用時間等。結(jié)果分析：包括R()的性質(zhì)、程序收斂情況、程序改進的方向等。每個同學必須獨立完成實驗（不能抄襲，否則兩人均為零分），實驗成績是該門課程成績的主要依據(jù)。六、實驗原理 1．R((S))的表示方法計算min和max是很容易的。；。當max時，R()=0。當minmax時，R()=min{I(U;V)：E(d)=}。在數(shù)學上，就是在約束條件：(1)(2)的約束下求平均信息量的條件極小值。為此引入待定常數(shù)S和u(u=1,2,...,r)，并作輔助函數(shù)(3)其中由得，為方便引入?yún)?shù)，則有(4)顯然(4)提供rs個方程，(2)提供r個方程，而(1)提供1各方程，共rs+r+1個方程；而有rs個未知數(shù)p(v|u)、r個未知數(shù)u及未知數(shù)S，共rs+r+1未知數(shù)，顯然可以求解。為方便起見，我們保留S作為參數(shù)。這樣得到：(5)(6)(7)(8)很容易得到，即S是率是失真函數(shù)的導數(shù)。當S-時，(S)min；參量S是的遞增函數(shù)，當從min到max逐漸增大時，S將隨增大而增大，當=max時,S達到最大值Smax<0。對Smax的求解較麻煩，必須解非線性方程。為了簡單我們不求Smax。如果r=s，即信源和信宿的符號集相同，則很容易通過(6)式求得u，進而通過(5)式求得p(v)。從而通過(7)(8)式劃出率失真函數(shù)曲線。2．R((S))的迭代計算 但一般情況下，rs，則只能通過(6)先求得p(v)，這是一個非常復雜的方程。下面介紹R((S))的迭代方法計算方法和公式。 首先假設(shè)p(v)固定，與信道傳遞概率p(v|u)無關(guān)，則求極值得： (9) 再假定p(v|u)不變，而把p(v)當成變量，則求極值得： (10) 具體算法為：選擇絕對值相當大的負數(shù)S1。選定起始傳遞概率p(1)(v|u)=1/rs。通過(10)式求得P(1)(v)，再通過(9)式求得p(2)(v|u)。如此重復直到與D(S1)(n+1)相差較??；并且。與R(S1)(n+1)相差較小再選擇較大的S2直到Smax逼近于零為止。這樣就可以畫出R()曲線.七。實驗步驟：1.建立熵函數(shù)求解matlab文件functionh=Hp(p)%熵函數(shù)計算，輸入是概率p（標量或一維矢量）%輸出是h,與p同維度h=zeros(1,length(p))index=find(p>0&p<1);P=p(index);h(index)=-P.*log(P)-(1-P).*log(1-P);end2.根據(jù)實驗原理，很容易知道，delta的最小值和最大植：delta_min=sum(Pu.*min(D',[],1));delta_max=min(Pu*D); 建立求率失真度的文件：%率失真函數(shù)計算function[R,delta]=R_delta(Pu,D,eps)%Pu信源概率矢量，D失真測度矩陣（rxs階）%信源Ur，信宿Vs[r,s]=size(D);%delta的最小值和最大植delta_min=sum(Pu.*min(D',[],1));delta_max=min(Pu*D);R=[];delta=[];%給定初始值%P(u,v)=P(V|U)P=ones(r,s)/s;%信道模型SS=100::-100;SS=-exp(SS);forS=SSPv=Pu*P;Ed0=sum(Pu*(P.*D));%Rs0=sum(Pu*(P.*log(P)))-sum(Pv.*log(Pv));Rs0=0;foru=1:rforv=1:sifP(u,v)~=0&Pu(u)~=0Rs0=Rs0+Pu(u)*P(u,v)*log(P(u,v)/Pv(v));endendendP=exp(S*D);fori=1:sP(:,i)=P(:,i)*Pv(i);endfori=1:rSumP=sum(P(i,:));P(i,:)=P(i,:)/SumP;endKm=50000;fork=1:KmPv=Pu*P;Edn=sum(Pu*(P.*D));%Rsn=sum(Pu*(P.*log(P)))-sum(Pv.*log(Pv));Rsn=0;foru=1:rforv=1:sifP(u,v)~=0&Pu(u)~=0Rsn=Rsn+Pu(u)*P(u,v)*log(P(u,v)/Pv(v));endendendP=exp(S*D);fori=1:sP(:,i)=P(:,i)*Pv(i);endfori=1:rSumP=sum(P(i,:));P(i,:)=P(i,:)/SumP;endifabs(Edn-Ed0)<eps&abs(Rsn-Rs0)<epsbreak;endEd0=Edn;Rs0=Rsn;endifk<KmR=[R,Rsn];delta=[delta,Edn];endend3.測試，以漢明失真測度矩陣為失真測度，分別繪制不同精度和不同輸入概率下的率失真函數(shù)曲線測試1：對于不同的精度測試程序如下clearallalpha=2;p=;D=[0,alpha;alpha,0];%漢明失真矩陣Pu=[p,1-p];%信源概率統(tǒng)計delta=[0:p/1000:p];R=Hp(p)-Hp(delta/alpha);%漢明失真，率失真函數(shù)理論值eps1=;%精度為eps1=[R1,delta1]=R_delta(Pu,D,eps1);eps2=;%精度為eps2=[R2,delta2]=R_delta(Pu,D,eps2);plot(delta,R,'b',delta1,R1,'r-.',delta2,R2,'m-*');legend('理論值','eps1=','eps2=')測試2：對于不同的概率輸入，測試程序如下：clearallalpha=2;p1=;p2=;D=[0,alpha;alpha,0];%漢明失真矩陣Pu1=[p1,1-p1];%信源概率統(tǒng)計Pu2=[p2,1-p2];%信源概率統(tǒng)計delta1=[0:p1/1000:p1];delta2=[0:p2/100