SDRookie
Junior Member level 2
Hi guys,
I am trying to simulate the 3rd intercept point with matlab. First, I generated 2 CW tones and it gave me the correct power of IM3. Then I generated a random noise signal and pass it through a bandpass filter, which gave me 2 tones with 1MHz bandwidth for each tone. The simulation results didn't give me the correct power as I expect. Can anyone help me please?
- - - Updated - - -
Modify the code make it more readable
I am trying to simulate the 3rd intercept point with matlab. First, I generated 2 CW tones and it gave me the correct power of IM3. Then I generated a random noise signal and pass it through a bandpass filter, which gave me 2 tones with 1MHz bandwidth for each tone. The simulation results didn't give me the correct power as I expect. Can anyone help me please?
Code:
close all;
clear;
%% Input Value
npoints=1e6;
rate=2e8;
R=50;
%% Bandlimit Gaussian Noise
freqbins = (0:npoints/2)/npoints*rate;
x=rand(1,npoints);
bpf1c=fdesign.bandpass('Fst1,Fp1,Fp2,Fst2,Ast1,Ap,Ast2',...
0.19,0.195,0.205,0.21,100,1,80);
bpf1=design(bpf1c,'equiripple');
bpf2c=fdesign.bandpass('Fst1,Fp1,Fp2,Fst2,Ast1,Ap,Ast2',...
0.29,0.295,0.305,0.31,80,1,80);
bpf2=design(bpf2c,'equiripple');
noise=filter(bpf1,x)+filter(bpf2,x);
%X=[1 exp(i*2*pi*rand(1,npoints/2-1)),1];
%X=fft(x,npoints);
% figure;
% pwelch(x,win,256,1024,rate,'power');
%% Zero out the frequency components outside the desired band
%X(find((freqbins < 19e6) | ((freqbins > 21e6)&(freqbins < 29e6))| (freqbins > 31e6))) = 0;
%NOISE= [X(1:npoints/2) conj(fliplr(X(1:npoints/2)))];
%NOISE=[X conj(fliplr(X(2:end-1)))];
%noise= ifft(NOISE);
% tc=gauspuls('cutoff',20e6,0.01,[],-40);
% t=-tc:1e-8:tc;
% noise=gauspuls(t,20e6,0.01);
[P,F]=periodogram(noise,[],npoints,rate,'power');
helperFrequencyAnalysisPlot2(F,10*log10(P)+30,'Frequency in Hz',...
'Power spectrum(dBm))',[],[],[0 1e8]);
%% non-linearity
G_db=20;
a1=(10^(G_db/10)*R)^0.5;
IP3_dBm=20;
IP3=10^(IP3_dBm/10)/1000;
a3=2*a1/(3*IP3*R);
IP3_expected=10*log10(2/(3*R)*1000)+10*log10(a1/a3);
%output=a1*noise+a3*noise.^3;
output_i=a3/a1*(noise.^3);
output_norm=output_i/((R)^0.5);
[Po,Fo]=periodogram(output_i,[],npoints,rate,'power');
helperFrequencyAnalysisPlot2(Fo,10*log10(Po)+30,'Frequency in Hz',...
'Power spectrum(dBm))',[],[],[0 1e8]);
%% Verify
IPj1=bandpower(noise,rate,[19e6 21e6]);
IPj1_dBm=10*log10(IPj1/R)+30;
IPj2=bandpower(noise,rate,[29e6 31e6]);
IPj2_dBm=10*log10(IPj2/R)+30;
IM3=bandpower(output_i,rate,[7e6 13e6]);
IM3_dBm=10*log10(IM3/R)+30;
IM3_dBm_expect=2*IPj1_dBm+IPj2_dBm-2*IP3_dBm;- - - Updated - - -
Modify the code make it more readable
Code:
close all;
clear;
%% Input Value
npoints=1e6;
rate=2e8;
R=50;
%% Bandlimit Gaussian Noise
freqbins = (0:npoints/2)/npoints*rate;
x=50*rand(1,npoints); %random noise
%bandpass filter
bpf1c=fdesign.bandpass('Fst1,Fp1,Fp2,Fst2,Ast1,Ap,Ast2',...
0.19,0.195,0.205,0.21,100,1,80);
bpf1=design(bpf1c,'equiripple');
bpf2c=fdesign.bandpass('Fst1,Fp1,Fp2,Fst2,Ast1,Ap,Ast2',...
0.29,0.295,0.305,0.31,80,1,80);
bpf2=design(bpf2c,'equiripple');
noise=filter(bpf1,x)+filter(bpf2,x); % noise with 2 tone
[P,F]=periodogram(noise,[],npoints,rate,'power');
helperFrequencyAnalysisPlot2(F,10*log10(P)+30,'Frequency in Hz',...
'Power spectrum(dBm))',[],[],[0 1e8]);
%% non-linearity
G_db=20; %power gain in db
a1=(10^(G_db/10)*R)^0.5; %linear voltage gain
IP3_dBm=20; %ip3
IP3=10^(IP3_dBm/10)/1000;
a3=2*a1/(3*IP3*R); %3rd order non linearity gain
IP3_expected=10*log10(2/(3*R)*1000)+10*log10(a1/a3);
%output=a1*noise+a3*noise.^3;
output_i=noise+a3/a1*(noise.^3); %refer back to input value
output_norm=output_i/((R)^0.5);
[Po,Fo]=periodogram(output_i,[],npoints,rate,'power');
helperFrequencyAnalysisPlot2(Fo,10*log10(Po)+30,'Frequency in Hz',...
'Power spectrum(dBm))',[],[],[0 1e8]);
%% Verify
%power of jammer1
IPj1=bandpower(noise,rate,[18e6 22e6]);
IPj1_dBm=10*log10(IPj1/R)+30;
%power of jammer2
IPj2=bandpower(noise,rate,[28e6 32e6]);
IPj2_dBm=10*log10(IPj2/R)+30;
%power of IM3
IM3=bandpower(output_i,rate,[7e6 13e6]);
IM3_dBm=10*log10(IM3/R)+30;
%expect IM3 value in theory
IM3_dBm_expect=2*IPj1_dBm+IPj2_dBm-2*IP3_dBm;