首先FM調制時屬于非線性調制中的一種，常見的非線性調制時FM(頻率調制)，PM(相位)調制

FM調制：就是利用調制信號的變化去改變載波的頻率，幅度保持不變，即載波信號的頻率隨我們調制信號的幅度變化而變化。

首先我們先介紹一些通信原理上的基本得概念：

現在FM調制就是我們將我們的基帶信號替換位其φt,但是我們會在其基帶信號的前面乘以一個調頻系數，也叫調頻靈敏度。這就通過我們基帶信號的振幅，改變其載波的頻率，這樣就是實現了頻率調制。下面我們假設基帶信號為一個余弦波：

上述就是其FM 調制的數學式子，這就實現了其FM 調制

這里說明上圖中的積分轉化為求和，因為積分是我們的連續域的，而我們的這都是數字信號，都是連續域，這就相當于我們的信號與系統中的傅里葉變化對應連續域，Z變化對應離散域。所以連續域的積分對應到離散域就是求和。

IQ調制：現代通信系統為了使頻譜利用率更高，用了許多種矢量調制,IQ調制中，兩個正交信號(頻率相同，相位相差90 °的載波，一般用Sin和Cos表示)與I(In-Phase，同相分量)、Q(Quadrature Phase，正交分量)兩路信號分別調制后一起發射，從而提高頻譜利用率。

所以接下來我們要將上式轉為IQ兩路，然后利用天線進行發送

我們I路信號為：cos ，圖中的cos(wt) 為載波信號

我們Q路信號為：sin , 圖中的sin(wt) 為載波信號

FM解調：

相干解調：相干解調也叫同步檢波，它適用于所有線性調制信號的解調。實現相干解調的關鍵是接收端要恢復出一個與調制載波嚴格同步的相干載波。相干解調是指利用乘法器，輸入一路與載頻相干(同頻同相)的參考信號與載頻相乘。

非相干解調：通信接收端從已調高頻信號中恢復出原始數字基帶信號時，采用的非相干解調方式，相對于相干解調方式，是指不需要提取載波信息的一種解調方法。

所謂相干，泛泛地說就是相互干擾;相干解調是指利用乘法器，輸入一路與載頻相干(同頻同相)的參考信號與載頻相乘。

這里的解調是運用的非相干解調：

我們將上述的調制信號分別乘以coswt,-sinwt,就可以得到相應的IQ路的數據，然后分別通過一個低通濾波器，濾除高頻分量，就可以的得到以基帶信號為頻率的I路Q路載波信號，然后兩個再相比，就可以得到其正切值，再求反正切就可以了的到其角度，角度就是我們的基帶信號。

這里的微分在離散域其實就是我們的相連的數據求差分。所以我們直接就可以求得基帶信號m.

FM調制matlab實現。

說明：因為后續會用硬件FPGA 實現，所以matlab 實現實現算法也會像盡可能的逼近其硬件。

思路：

產生基本的信號：基帶信號m(n),載波cos(wt),sin(wt)的實現

我們這里實現信號發生器是用的查表法，即FPGA，MATLAB 產生DDS,這里就不說明其方法了，相關信息可以看我發的其他博客，里面有詳細的教程，以及原理介紹。

下面就是FM調制解調代碼：

clc;

clear all;

%*****************************************

fs = 16e6; %采樣率 載波的采樣率

fc = 1e6; %載波中心頻率

df = 75e3; %最大頻偏

fm = 16e3; %音頻的采樣率

kf = (df*2^32/fs)/32767;

AC= 1024;%幅度

%由于我們的載波的采樣率是我們的音頻的采樣率的1000倍，但是我們在進行FM調制時

%會用到兩個相乘，所以相當于一個語音信號的點，要被采樣1000次，并且要求數據長度一致

%所以我們可以將語音信號的每一個樣點復制1000次

%read pcm file

%*********************

%read pcm file

file= fopen('v3edu.pcm','r');

m=fread(file,'int16');

m_len=length(m)*fs/fm;%按照載波采樣率的長度計算

figure(1);

plot(m);title('原始的語音信號');

%*****************************

%generate 1Mhz

w=fc*2^32/fs; %頻率控制字，32位的相位累加器， f_out = fs*w/2^N :

% 載波的ROM

n=0:1/1024:1023/1024;

s_rom=sin(2*pi*n);

c_rom=cos(2*pi*n);

%plot(s_rom);

%接下來產生載波地址，

w_r=0;%相位累加器

rrom_addr=0;%

cw_sin=zeros(1,m_len);%這個相當于是我們的查找表，一個周期的正弦波

cw_cos=zeros(1,m_len);

for i=1:m_len

w_r = w_r + w;

if(w_r > 2^32) % 做32位累加器的溢出判斷

w_r= w_r - 2^32;

end

rrom_addr=round(w_r/2^22);%讀查找表的地址

if rrom_addr == 0%matlab計數時從1開始計數

rrom_addr =1;

end

cw_sin(i)=s_rom(rrom_addr);%%%%%%載波完成

cw_cos(i)=c_rom(rrom_addr);%載波

end

figure(2);

subplot(2,1,1);plot(cw_sin(1:3000),'r');title('sin載波信號');

subplot(2,1,2);plot(cw_cos(1:3000),'b');title('cos載波信號');

%給音頻數據升采樣，與載波長度一致，后面要相乘，所以要保證信號的長度是一樣的

m_t=zeros(1,m_len);

for i=1:length(m)

for j=1:fix(fs/fm)

m_t((i-1)*fix(fs/fm)+j)=m(i);%同一個16k的采樣點復制1000次這樣就是16M采樣點

end

end

%*************************

%調制 sin 和cos 的信號

w_r=0;%相位累加器

rrom_addr=0;%

dac_i=zeros(1,m_len);%cos，I路信號

dac_q=zeros(1,m_len);%sin Q路信號

for i=1:m_len

w_r = w_r + kf*m_t(i);

if(w_r > 2^32) % 做32位累加器的溢出判斷

w_r= w_r - 2^32;

elseif(w_r <0)

w_r = w_r + 2^32; % 負的溢出時

end

rrom_addr=round(w_r/2^22);%讀查找表的地址

if rrom_addr == 0

rrom_addr =1;

end

dac_q(i)=AC*s_rom(rrom_addr);

dac_i(i)=AC*c_rom(rrom_addr);%載波

end

%****************************

%正交調制

s_t=zeros(1,m_len);

for i=1:m_len

s_t(i) = dac_i(i)*cw_cos(i) + dac_q(i)*cw_sin(i)*(-1);

end

figure(3);

plot(s_t);title('調制完成信號');

調制結束

MATLAB 濾波器的實現：

我們實現一個簡單的低通濾波器，只需要將我們我們需要的濾波器參數輸入Filter Solutions中，然后復制器抽頭系數，然后與我們的輸入數據相卷積就可以了，下面的代碼有注釋.

%%%%%%%%%%%%正交解調%%%%%%%%%%%%%%%%

i_data = zeros(1,m_len);

q_data = zeros(1,m_len);

for i = 1 : m_len

i_data(i) = s_t(i) * cw_cos(i);

q_data(i) = (-1)* s_t(i) * cw_sin(i);

end

%%%%%%%%%%%%%%進行低通濾波%%%%%%%%%%%%%%%%

%濾波器抽抽頭系數，

NUM = [-7.894e-05, -2.483e-04, -4.516e-04, -7.089e-04, -1.035e-03, -1.435e-03, -1.903e-03, -2.417e-03, -2.938e-03, -3.415e-03, -3.778e-03, -3.948e-03, -3.837e-03, -3.357e-03, -2.42e-03, -9.524e-04, 1.106e-03, 3.791e-03, 7.113e-03, 1.105e-02, 1.554e-02, 2.051e-02, 2.583e-02, 3.135e-02, 3.691e-02, 4.233e-02, 4.742e-02, 5.199e-02, 5.588e-02, 5.894e-02, 6.105e-02, 6.213e-02, 6.213e-02, 6.105e-02, 5.894e-02, 5.588e-02, 5.199e-02, 4.742e-02, 4.233e-02, 3.691e-02, 3.135e-02, 2.583e-02, 2.051e-02, 1.554e-02, 1.105e-02, 7.113e-03, 3.791e-03, 1.106e-03, -9.524e-04, -2.42e-03, -3.357e-03, -3.837e-03, -3.948e-03, -3.778e-03, -3.415e-03, -2.938e-03, -2.417e-03, -1.903e-03, -1.435e-03, -1.035e-03, -7.089e-04, -4.516e-04, -2.483e-04, -7.894e-05];

adc_i = conv(i_data,NUM); %低通濾波

adc_q = conv(q_data,NUM);%%%%

%%%%%%%%%%%%%FM解調%%%%%%%%%%%%%%%%%%%%%%

%重新檢測數據長度

c_len = length(adc_i);

cr = zeros(1,c_len);

cj = zeros(1,c_len);

for i = 2: c_len %必須從2開始，因為有I-1

cr(i) = adc_i(i)*adc_i(i) +adc_q(i)*adc_q(i);

cj(i) = adc_i(i-1) * adc_q(i) - adc_i(i)* adc_q(i-1);%I(n-1)*Q(n) -I(n)*Q(n-1)

end

angle = zeros(1,c_len);

for i = 1:c_len

if cr(i) == 0

angle(i) = 0;

else

angle(i) = (cj(i)/cr(i))*2 ^16; %FPGA 中沒有小數，要進行量化

end

end

%%%%%%%%%%%%%降采樣%%%%%%%%%%%%%%%%%%%%%%%%%

%降采樣，以均值方式，16M，降到16k

d_len = length(angle);

sum = 0;

cnt = 0;

demout = zeros(1,fix(d_len*(fm/fs)));

d_cnt = 1;

for i= 1:d_len

sum = sum + angle(i);

if cnt == fs/fm %完成1000次的累加

demout(d_cnt) = sum /cnt;

d_cnt = d_cnt + 1;

cnt = 0;

sum= 0;

end

cnt = cnt + 1;

end

figure(4);

plot(demout);title('解調后的信號');

最后結果圖：