断更了很久，一方面是因为之后开始准备期末考试和课程设计，另一方面也是因为自己因为放假堕落了（好吧，这个是主要原因）。
课程设计做了挺久，整合了很多，参考了不少内容才做出来。服务学弟学妹了，就把我的这个发出来好了。
其中，信源为自己构造的时间函数，数字化方式为PCM编码，基带码为CMI码，信道码为汉明码，调制方式为ASK调制与解调，信道类型为衰落信道，接收端与发送端对应。
1.main函数
%empty the processing
clc;
close all;
clear all;%original signal
clear;
fs=8000;%抽样频率
fc=100;
dt=1/fs;
Ts=dt;
t=0:dt:0.02-dt;
y=0.3*cos(2*pi*fc*t)+0.5*sin(4*pi*fc*t)+0.2*sin(8*pi*fc*t);%输入模拟信号
figure(1);
plot(t,y);
ylim([-2,2]);
title('源信号');%PCM modulation
n=length(t);
y_PCM=PCMENC2(y,n);
figure(2);
stairs(y_PCM);
%zoom xon;
xlim([0 100]);
ylim([-0.2,1.2]);
title('PCM编码信号');%CMI modulation
y_CMI = CMI_coding(y_PCM) ;
t2 = 0:Ts:41*Ts;
figure(3)
stairs(y_CMI);
xlim([0 200]);
ylim([-0.2 1.2]);
%zoom xon;axis([0 length(y_CMI) -0.2 1.2]);
title('CMI编码信号');%Hamming modulation
y_Hamming=encode(y_CMI,7,4,'hamming') ;
figure(4)
stairs(y_Hamming);
xlim([0 200]);
ylim([-0.2 1.2]);
%zoom xon;axis([0 length(y_Hamming) -0.2 1.2]);
title('Hamming编码');%ASK modulaiton and rayleigh channel
y_DASK = ASK(y_Hamming);%Hamming demodulation
y_DHamming=decode(y_DASK,7,4,'hamming');
figure(8)
stairs(y_DHamming);
xlim([0 200]);
ylim([-0.2 1.2]);
%zoom xon;axis([0 length(y_DHamming) -0.2 1.2]);
title('Hammming译码输出信号（CMI）');%CMI demodulation
y_DCMI = CMI_decoding(y_DHamming);
figure(9)
stairs(y_DCMI);
xlim([0 100]);
ylim([-0.2 1.2]);
%zoom xon;axis([0 length(y_DCMI) -0.2 1.2]);
title('CMI译码输出信号（PCM信号）');%PCM demodulation
y_DPCM=PCMDEC2(y_DCMI,n);
figure(10);
plot(t,y_DPCM);
ylim([-2 2]);
title('解码输出');%original signal modulation
figure(11)
plot(t,y);ylim([-2 2]);hold on;
plot(t,y_DPCM);ylim([-2 2]);
title('发送和接收信号对比，蓝色为源信号，红色为接收信号');%signal distortion
di=0;
for k=1:nd=(y(k)-y_DPCM(k))^2/n;di=di+d;
end
di=sqrt(di);
disp('失真度为 ');
disp(di);2.PCM编码函数
function a=PCMENC2(y,n) %S为输入信号
a=zeros(1,8*n);
for i=1:na(8*(i-1)+1)=(sign(y(i))+1)/2;s=round(abs(y(i))*2047);if s<16a(8*(i-1)+2:8*(i-1)+4)=[0 0 0];temp=0;del=1;elseif s<32a(8*(i-1)+2:8*(i-1)+4)=[0 0 1];temp=16;del=1;elseif s<64a(8*(i-1)+2:8*(i-1)+4)=[0 1 0];temp=32;del=2;elseif s<128a(8*(i-1)+2:8*(i-1)+4)=[0 1 1];temp=64;del=4;elseif s<256a(8*(i-1)+2:8*(i-1)+4)=[1 0 0];temp=128;del=8;elseif s<512a(8*(i-1)+2:8*(i-1)+4)=[1 0 1];temp=256;del=16;elseif s<1024a(8*(i-1)+2:8*(i-1)+4)=[1 1 0];temp=512;del=32;elseif s<2048a(8*(i-1)+2:8*(i-1)+4)=[1 1 1];temp=1024;del=64;endb=zeros(1,4);m=4;s1=fix((s-temp)/del);while s1>0b(1,m)=mod(s1,2);s1=(s1-b(1,m))/2;m=m-1;enda(8*(i-1)+5:8*i)=b;
end3.CMI编码函数
function y = CMI_coding(y1)a=reshape(y1',1,numel(y1));
n=0;
m=1;
b=zeros(1,length(a));
c=zeros(1,length(a));
for k=1:length(a)m=mod(n,2); %CMI编码结果“00”和“11”交替的判决语句;if (a(k)==0) %如果输入码为“0”，则进行下一步;b(k)=0; %CMI编码输出的前一位为“0”;c(k)=1; %CMI编码输出的后一位为“1”;endif (a(k)==0)continue;endif (a(k)==1 && m==1) %如果输入码为“1”且判决码为1，则进行下一步;b(k)=1; %CMI编码输出的前一位为“1”;c(k)=1; %CMI编码输出的后一位为“1”;n=n+1; %判决辅助码加“1”endif (a(k)==1 && m==0) %如果输入码为“1”且判决码为0，则进行下一步;b(k)=0; %CMI编码输出的前一位为“0”;c(k)=0; %CMI编码输出的后一位为“0”;n=n+1; %判决辅助码加“1”endend%display('Display encode result:');
y=zeros(1,2*length(a));
for k=1:length(a)if (b(k)==1 &&c(k)==1) %如果编码后的前一位和后一位都为“1”，则进行下一步;y(1,2*k-1)=1;y(1,2*k)=1;%fprintf(' 1 1 ') ; %CMI编码输出为"11";else if (b(k)==0 &&c(k)==0) %如果编码后的前一位和后一位都为“0”，则进行下一步;y(1,2*k-1)=0;y(1,2*k)=0; %fprintf(' 0 0 ') ; %CMI编码输出为"00";else y(1,2*k-1)=0;y(1,2*k)=1; % fprintf(' 0 1 ') ; %CMI编码输出为"01";endend
end4.ASK编译码函数及衰落信道
function y_DASK =ASK(code)codn=length(code);
fc=10;
fs=fc*6;
code_len=8;
for i=1:codn % 产生数字基带信号x((i-1)*code_len+1:code_len*i)=code(i);
end
car=cos(2*pi*fc/fs*(0:length(x)-1));
y_ASK=x.*car;
figure(5);
plot(y_ASK);
xlim([0 800]);
ylim([-5 5]);
%zoom xon;axis([0 length(y_ASK) -1.2 1.2]);
title('ASK调制信号');%衰落信道
%y_rayleigh=Rayleigh(y_ASK);
%figure(6);
%plot(y_rayleigh);zoom xon;
%axis([0 800 -2 2]);
%axis([0 length(y_rayleigh) -2 2]);
%title('衰落信道');
chan=rayleighchan(1/15000,0);
%1/15000代表我们的采样频率，也就是我一个点代表的是多长一段时间单位。
%0就是最大多普勒频移
y_rayleigh=filter(chan,y_ASK);
figure(6);
plot(abs(y_rayleigh));zoom xon;
axis([0 800 -5 5]);
%axis([0 length(y_rayleigh) -5 5]);
title('衰落信道'); %解调
yn3=y_rayleigh.*car;
for i=1:length(yn3)if(yn3(i)>=0.02)yn3(i)=1;else yn3(i)=0;end
end
y_DASK=zeros(1,codn);
i1=1;
for i=1:8:length(yn3)y_DASK(i1)=yn3(i);i1=i1+1;
end
figure(7);
%plot(y_DASK);
stairs(y_DASK);
xlim([0 200]);
ylim([-0.2 1.2]);
%zoom xon;axis([0 length(y_DASK) -0.2 1.2]);
title('ASK解调信号（Hamming）');5.CMI译码函数
function y_DCMI=CMI_decoding(y)a=reshape(y',1,numel(y));
b=zeros(1,length(a)/2);
for i=1:length(a)/2 b(i)=a(i*2-1);
end
c=zeros(1,length(a)/2);
for i=1:length(a)/2 c(i)=a(i*2);
end
y_DCMI=zeros(1,length(a)/2);
%display('Display decode result:');
for k=1:length(a)/2if(b(k)==1 &&c(k)==1) %如果CMI码的前序列和后序列为"1",则运行下一步;y_DCMI(k)=1;%fprintf(' 1 '); %输出显示"1";endif (b(k)==0 &&c(k)==0) %如果CMI码的前序列和后序列都为"0",则运行下一步;y_DCMI(k)=1;%fprintf(' 1 '); %输出显示"1";endif (b(k)==0 &&c(k)==1) %如果CMI码的前序列为"0"且后序列为"1",则运行下一步;y_DCMI(k)=0;%fprintf(' 0 '); %输出显示"0";end
end6.PCM译码函数
function b=PCMDEC2(code,n)for i=1:nif code(8*(i-1)+2:8*(i-1)+4)==[0 0 0]temp=0;del=1;elseif code(8*(i-1)+2:8*(i-1)+4)==[0 0 1]temp=16;del=1;elseif code(8*(i-1)+2:8*(i-1)+4)==[0 1 0]temp=32;del=2;elseif code(8*(i-1)+2:8*(i-1)+4)==[0 1 1]temp=64;del=4;elseif code(8*(i-1)+2:8*(i-1)+4)==[1 0 0]temp=128;del=8;elseif code(8*(i-1)+2:8*(i-1)+4)==[1 0 1]temp=256;del=16;elseif code(8*(i-1)+2:8*(i-1)+4)==[1 1 0]temp=512;del=32;elseif code(8*(i-1)+2:8*(i-1)+4)==[1 1 1]temp=1024;del=64;endb1=0;for m=5:8b1=2*b1+code(8*(i-1)+m);endb(i)=(2*code(8*(i-1)+1)-1)*round(temp+(b1+0.5)*del)/2047;
end
