Files
imdd_silas/projects/FSO_transmission/first_analysis.m
silas (home) dd6094d51e chirped WDM analsyis
minor stuff in FSO link with Magnus
2026-01-12 08:12:42 +01:00

172 lines
6.9 KiB
Matlab
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
base = "C:\Users\Silas\Nextcloud\Dokumente\02_Ablage_Office\FSO_FP_QCL_60umUTC";
mode = 0; %0 oder 1
M = 2;
all_files = dir(fullfile(base, "**/*.mat"));
if M == 2
tx_data = load("C:\Users\Silas\Nextcloud\Dokumente\02_Ablage_Office\FSO_FP_QCL_60umUTC\14G_PAM2\tx_info\tx_info_PAM2_14Gbd0.75RRC.mat");
filename = fullfile(base, "14G_PAM2\M=2_Rs=1.4e10_Fs=8e10_I=265mA_RoP=46.3mW_L=31m_PS=RRC_rolloff=0.75_Mode=Rise.mat");
elseif M == 4
tx_data = load("C:\Users\Silas\Nextcloud\Dokumente\02_Ablage_Office\FSO_FP_QCL_60umUTC\6G_PAM4\tx_info\tx_info_PAM4_6Gbd0.6RRC.mat");
filename = fullfile(base, "6G_PAM4\M=4_Rs=6e9_Fs=8e10_I=255mA_RoP=42.3mW_L=31m_PS=RRC_rolloff=0.6_Mode=Rise.mat");
end
if mode == 1
[f, p] = uigetfile(fullfile(base, "**/*.mat"));
if f~=0
filename = fullfile(p,f);
end
end
datas = load(filename);
%%
str = filename;
M_ = str2double(regexp(str, 'M=([^_]+)', 'tokens', 'once'));
assert(M==M_);
fsym = str2double(regexp(str, 'Rs=([^_]+)', 'tokens', 'once'));
fs = str2double(regexp(str, 'Fs=([^_]+)', 'tokens', 'once'));
I = sscanf(char(regexp(str, 'I=([^_]+)', 'tokens', 'once')), '%f');
rop = sscanf(char(regexp(str, 'RoP=([^_]+)', 'tokens', 'once')), '%f');
L = sscanf(char(regexp(str, 'L=([^_]+)', 'tokens', 'once')), '%f');
pulseshape = string( regexp(str, 'PS=([^_]+)', 'tokens', 'once'));
rolloff = str2double(regexp(str, 'rolloff=([^_]+)', 'tokens', 'once'));
mode = string( regexp(str, 'Mode=([^\.]+)', 'tokens', 'once'));
%%
% Tx data
Bits = Informationsignal(tx_data.tx_data,"fs",fsym);
Symbols = Informationsignal(real(tx_data.tx_PAM_sym),"fs",fsym);
mapping_style = M==4; % Pam2 is like move-it; PAM-4 is different, same mapping like ETH peopled used in Zurich... hence the "eth_style" argument here and there
PM = PAMmapper(M,0,"eth_style",mapping_style); % one should rename "eth style" as this is simply a different mapping scheme
Symbols_ = PM.map(Bits) .* PM.scaling;
assert(isequal(Symbols.signal,Symbols_.signal));
Bits_ = PM.demap(Symbols);
[bits,errors,ber,errorIndice] = calc_ber(Bits_.signal,Bits.signal);
assert(ber == 0);
%% For comparison, apply pulsef on Tx Symbols
Pform = Pulseformer("fsym",fsym,"fdac",fs,"pulse","rc","pulselength",16,"alpha",rolloff);
Digi_sig_compare = Pform.process(Symbols);
MF = Pulseformer("fsym",fsym,"fdac",2*fsym,"pulse","rrc","pulselength",16,"alpha",rolloff);
Rx_sig_compare = MF.process(Digi_sig_compare);
%%
% Rx Data
traceData = datas.tr.lastData(2).trace.ch3;
%FYI: Voltage=(RawDataYReference)×YIncrement+YOrigin
scoperead_volts = (traceData.RawData - traceData.YReference) * traceData.YIncrement + traceData.YOrigin;
demystified = isequal(traceData.YData,scoperead_volts);
assert(demystified);
Scope_sig = Electricalsignal(traceData.YData,"fs",fs);
Scope_sig.plot("displayname",'raw','fignum',100);
Scope_sig.spectrum("displayname",'raw','fignum',101)
% 1) matched filter
% pulse is symmetric, hence we can use pulsef firectly as matched filter.
% It feels off (bit I think correct) that the fsym is now the output freq.!!
% -> output 2 sps to omit timing recovery!?
Pform = Pulseformer("fsym",fsym,"fdac",2*fsym,"pulse","rrc","pulselength",16,"alpha",rolloff,"matched",1);
Rx_matched = Pform.process(Scope_sig);
Rx_matched.spectrum("displayname",'Signal after matched filter','fignum',1);
%%
sys = comm.SymbolSynchronizer('TimingErrorDetector', 'Gardner (non-data-aided)', ...
'SamplesPerSymbol', 2, ...
'DampingFactor', 0.7, ...
'NormalizedLoopBandwidth', 0.01);
Rx_symbolsync = Rx_matched;
[Rx_symbolsync.signal, timing_error] = sys(Rx_matched.signal);
plot(timing_error); % If this is a ramp, you have drift!
%% timing sync -> at this point we still have no symbol timing recovery, we
% % try to do this with 2sps EQ!
[~,Rx_synced_cell,inverted,sequenceFound,sequenceStarts] = Rx_symbolsync.tsynch("reference", Symbols, "fs_ref", fsym, "debug_plots", 1);
%% not working..
Rx_synced = Rx_synced_cell{1};
len_tr = 4096*2;
mu_ffe1 = 0.0001;
mu_ffe2 = 0.0008;
mu_ffe3 = 0.001;
mu_dc = 0.005;
mu_ffe = [mu_ffe1 mu_ffe3 mu_ffe3];
mu_dfe = 0.0004;
duob_mode = db_mode.no_db;
Rx_synced.plot("displayname",'RX: Matched+Sync+2sps','fignum',2);
Digi_sig_compare.normalize("mode","rms").spectrum("displayname",'Tx: RC-shaped','fignum',1,'normalizeTo0dB',0);
Rx_sig_compare.normalize("mode","rms").spectrum("displayname",'Tx: RC-shaped + matched filtered ','fignum',1,'normalizeTo0dB',0);
Rx_synced.normalize("mode","rms").spectrum("displayname",'RX: matched filtered + synced','fignum',1,'normalizeTo0dB',0);
if M == 2
ber_in_paper = 10^(-2.6); %fig 3a) 4 Gb/s MWIR FSO Transmission using Directly Modulated QCL and an Uncooled UTC-PD at Room-Temperature
elseif M == 4
ber_in_paper = 10^(-2.5);
end
%% -------------------- FFE --------------------
% requires some more digging what is going on :-)
eq_ffe = EQ("Ne",[50, 5, 5],"Nb",[2,0,0], ...
"training_length",len_tr,"training_loops",5,"dd_loops",5, ...
"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005, ...
"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
ffe_results = ffe(eq_ffe,M,Rx_synced,Symbols,Bits, ...
"precode_mode",duob_mode,'showAnalysis',0,"postFFE",[], ...
"eth_style_symbol_mapping",mapping_style);
% ffe_results.metrics.print
fprintf('My EQ: %.1e \n',ffe_results.metrics.BER);
fprintf('Paper: %.1e \n \n',ber_in_paper);
%% -------------------- VNLE + MLSE --------------------
pf_ncoeffs = 1;
eq_v = EQ("Ne",[100, 5, 5],"Nb",[0, 0, 0], ...
"training_length",len_tr,"training_loops",5,"dd_loops",5, ...
"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005, ...
"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
pf_ = Postfilter("ncoeff",pf_ncoeffs,"useBurg",1);
mlse_ = MLSE("duobinary_output",0,'M',M,'trellis_states',PAMmapper(M,0,"eth_style",mapping_style).levels);
[vnle_results, mlse_results] = vnle_postfilter_mlse(eq_v, pf_, mlse_, M, Rx_synced, Symbols, Bits, ...
"precode_mode", duob_mode, 'showAnalysis', 1, "postFFE", [], "eth_style_symbol_mapping", mapping_style);
mlse_results.metrics.print("description",'MLSE')
fprintf('My EQ: %.1e \n',mlse_results.metrics.BER);
fprintf('Paper: %.1e \n \n',ber_in_paper);
%% -------------------- DB target --------------------
mlse_db_ = MLSE("DIR",[1,1],"duobinary_output",0,"M",M,'trellis_states',PAMmapper(M,0).levels);
eq_ = EQ("Ne",[50, 5, 5],"Nb",[0,0,0],"training_length",len_tr,"training_loops",5,"dd_loops",5, ...
"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
dbt_results = duobinary_target(eq_,mlse_db_, M, Rx_synced, Symbols, Bits, ...
"precode_mode", duob_mode, 'showAnalysis', 0, "postFFE", [],"eth_style_symbol_mapping",mapping_style);
dbt_results.metrics.print("description",'Duobinary');
mlse_results.metrics.print
fprintf('My EQ: %.1e \n',dbt_results.metrics.BER);
fprintf('Paper: %.1e \n \n',ber_in_paper);