diff --git a/Classes/00_signals/Electricalsignal.m b/Classes/00_signals/Electricalsignal.m index 2d28cfe..95beb14 100644 --- a/Classes/00_signals/Electricalsignal.m +++ b/Classes/00_signals/Electricalsignal.m @@ -28,8 +28,10 @@ classdef Electricalsignal < Signal function pow = power(obj) - pow = mean( abs(obj.signal).^2 ) ; - pow = pow2db(pow)+30; %dbm + % Power of an electrical signal + R = 50; + pow = mean(abs(obj.signal).^2) / R; % Power in watts = V^2 / R + pow = 10*log10(pow) + 30; % Power in dB +30 = dBm end @@ -47,6 +49,66 @@ classdef Electricalsignal < Signal end + function obj = normalize(obj,options) + + arguments + obj Electricalsignal + options.mode normalization_mode = normalization_mode.rms + end + + switch options.mode + + case normalization_mode.milliwatt + + curpow = sqrt(mean(obj.signal.^2)/50); %watt + tgtpow = sqrt(1e-3); %watt -> wir wollen milliwatt + + scling = tgtpow/curpow; + + obj.signal = obj.signal*scling; + + case normalization_mode.rms + + obj.signal = obj.signal / sqrt(mean(obj.signal.^2)); + + case normalization_mode.oneone + + obj.signal = obj.signal - min(obj.signal); + obj.signal = obj.signal/max(abs(obj.signal)); + obj.signal = (2*obj.signal) - 1; + + end + + end + + function obj = setPower(obj,outputpower, mode) + + arguments + obj Electricalsignal + outputpower + mode power_notation + end + + switch mode + case power_notation.W + + case power_notation.mW + outputpower_w = outputpower/1000; %mW to Watt + case power_notation.dBm + outputpower_w = 10^((outputpower-30)/10); % dBm to Watt + case power_notation.dBW + outputpower_w = 10^((outputpower)/10); % dBm to Watt + end + + curpow = sqrt(mean(obj.signal.^2)/50); %watt mit 50 ohm ref + tgtpow = sqrt(outputpower_w); %watt -> wir wollen milliwatt + + scling = tgtpow/curpow; + + obj.signal = obj.signal*scling; + + end + end diff --git a/Classes/00_signals/Informationsignal.m b/Classes/00_signals/Informationsignal.m index afa7462..03ccc00 100644 --- a/Classes/00_signals/Informationsignal.m +++ b/Classes/00_signals/Informationsignal.m @@ -28,7 +28,7 @@ classdef Informationsignal < Signal function pow = power(obj) - pow = mean(abs(obj.signal),"all") ; + pow = mean(abs(obj.signal.^2),"all") ; end diff --git a/Classes/00_signals/Signal.m b/Classes/00_signals/Signal.m index 8080d4c..758bb96 100644 --- a/Classes/00_signals/Signal.m +++ b/Classes/00_signals/Signal.m @@ -373,14 +373,14 @@ classdef Signal fsig = obj.fs; q = fsig/fsym; - if q > 10 - sig = abs(obj.signal).^2; + if q > 10 && isinteger(q) + sig = (obj.signal); else - sig = abs(obj.resample("fs_in",fsig,"fs_out",fsym*10).signal); + sig = (obj.resample("fs_in",fsig,"fs_out",fsym*30).signal); q = 10; end - figure(12) + figure() clf cursor = 200*q; @@ -390,6 +390,8 @@ classdef Signal cursor=cursor+q; s=s+1; end + + ylim([-3 3]); end diff --git a/Classes/01_transmit/AWG.m b/Classes/01_transmit/AWG.m index 47ecd23..9f2ba69 100644 --- a/Classes/01_transmit/AWG.m +++ b/Classes/01_transmit/AWG.m @@ -5,9 +5,10 @@ classdef AWG properties(Access=public) kover %oversampling factor e.g. 16 + upsampling_method repetitions %repeat the signal to generate a longer sequence? fdac %needed - normalize %want to normalize at first? either 0 or 1 + normalize2dac %want to normalize at first? either 0 or 1 bit_resolution %bit res. of quantizer (e.g. 5 bit) dac_min dac_max @@ -32,8 +33,9 @@ classdef AWG arguments options.kover = 16; + options.upsampling_method upsampling_mode options.repetitions = 1; - options.normalize = 1; + options.normalize2dac = 1; options.fdac = 92e9; options.bit_resolution = 5.5 options.dac_min = -0.5; @@ -44,6 +46,7 @@ classdef AWG options.lpf_type = 0; options.f_cutoff = 32e9; options.H_lpf Filter + end fn = fieldnames(options); @@ -59,9 +62,20 @@ classdef AWG len_in = length(signalclass_in.signal); + if signalclass_in.fs ~= obj.fdac + signalclass_in = signalclass_in.resample("fs_in",signalclass_in.fs,"fs_out",obj.fdac); + end + % 1-3. actual processing of the signal (normalize->quantize->sample hold) signalclass_in.signal = obj.process_(signalclass_in.signal); + + % cast the inform. signal to electrical signal + signalclass_in = Electricalsignal(signalclass_in,"fs",obj.fdac*obj.kover,"logbook",signalclass_in.logbook); + % normalize to 0dBm before applying the lowpass + %signalclass_in = signalclass_in.normalize("mode","milliwatt"); + signalclass_in = signalclass_in.setPower(12,"dBm"); + % 4. Apply LPF on the signal if obj.lpf_active if isa(obj.H_lpf,'Filter') @@ -74,10 +88,9 @@ classdef AWG signalclass_in = lpf.process(signalclass_in); end end + - % cast the inform. signal to electrical signal - signalclass_in = Electricalsignal(signalclass_in,"fs",obj.fdac*obj.kover,"logbook",signalclass_in.logbook); - + % append to logbook current_class = class(obj); lbdesc = ['AWG ', current_class , '// k_over:',num2str(obj.kover),'. f_dac:',num2str(obj.fdac*1e-9),'GHz. Resolution:',num2str(obj.bit_resolution),' bits.']; @@ -108,16 +121,17 @@ classdef AWG obj.signal_length = length(data_in); - if obj.normalize - % 0a Normalize the signal to 1 Vpp and set the amplitude of the signal + if obj.normalize2dac + % 0a Normalize the signal to full scale DAC range + data_in = data_in - min(data_in); data_in = data_in/(max(data_in)-min(data_in)); - else - % 0b Cut the Signal at -1 and 1 and scale to amplitude - data_in(data_in > 1) = 1; - data_in(data_in < -1) = -1; + data_in = data_in * (obj.dac_max-obj.dac_min); + data_in = data_in + obj.dac_min; end - % 1. Quantize the signal + % 1. Quantize the signal - Full Scale is between obj.dac_min + % and dac_max. If signal is smaller in between, you won't use + % the full bit-resolution. if obj.bit_resolution>0 elec_out = obj.quantization(data_in) ; else @@ -125,8 +139,16 @@ classdef AWG end % 2. Sample and hold + repeat (data_out: 1xsignal length) - elec_out = repmat(elec_out,obj.repetitions,obj.kover); - elec_out = reshape(elec_out',[],1); + if obj.upsampling_method == 1 + % just use matlab function + elec_out = resample(elec_out,obj.kover,1); + elseif obj.upsampling_method == 2 + % sample and hold + elec_out = repmat(elec_out,obj.repetitions,obj.kover); + elec_out = reshape(elec_out',[],1); + else + error('chosen upsampling method not implemented?'); + end % 3. Add skew (not implemented so far) if obj.skew_active @@ -145,7 +167,7 @@ classdef AWG if isreal(x_in) % shift signal and clip to quantizer intervall - x_in = min(max(x_in-obj.dac_min,0),obj.dac_max-obj.dac_min); + x_in = min( max(x_in-obj.dac_min,0), obj.dac_max-obj.dac_min); % quantize signal x_in = round((steps-1)/(obj.dac_max-obj.dac_min)*x_in); diff --git a/Classes/01_transmit/PAMmapper.m b/Classes/01_transmit/PAMmapper.m index 6ff22ec..ed015fe 100644 --- a/Classes/01_transmit/PAMmapper.m +++ b/Classes/01_transmit/PAMmapper.m @@ -24,6 +24,7 @@ classdef PAMmapper if isa(signal_in,'Signal') signal_in.signal = obj.map_(signal_in.signal); + signal_in = signal_in.normalize("mode","rms"); signal_in = signal_in.logbookentry(); out = signal_in; else @@ -57,7 +58,7 @@ classdef PAMmapper pam_sig=2*pam_sig-3; end - pam_sig = pam_sig .* 1/sqrt(5); + pam_sig = pam_sig/sqrt(5); case 6 @@ -85,7 +86,7 @@ classdef PAMmapper pam_sig=2*pam_sig-7; end - + pam_sig = pam_sig/sqrt(21); case 16 % 16-ASK: x1 = bitpattern(:,1); @@ -126,6 +127,7 @@ classdef PAMmapper elseif obj.unipolar==1 thres=[0.5,1.5,2.5]; end + thres = thres .* 1/sqrt(5); case 6 %PAM 6 @@ -139,6 +141,8 @@ classdef PAMmapper elseif obj.unipolar==1 thres=0.5:6.5; end + + thres=thres./sqrt(21); case 16 % 16-ASK diff --git a/Classes/01_transmit/Pulseformer.m b/Classes/01_transmit/Pulseformer.m index e51c601..a3fa2cd 100644 --- a/Classes/01_transmit/Pulseformer.m +++ b/Classes/01_transmit/Pulseformer.m @@ -88,7 +88,7 @@ classdef Pulseformer %Bau das Filter (hier rrc) racos_len = obj.pulselength*2; alpha = obj.rrcalpha; - h = rcosdesign(alpha,racos_len,sps); + h = rcosdesign(alpha,racos_len,sps,"normal"); h = h./ max(h); end @@ -118,8 +118,8 @@ classdef Pulseformer % %Apply Filter using Matlab build in fctn. - h = rcosdesign(alpha,racos_len,sps); - h = h./ max(h); +% h = rcosdesign(alpha,racos_len,sps); +% h = h./ max(h); %data_out_ = upfirdn(data_in,h,up,dn); % % %cut signal, which is longer due to fir filter @@ -128,8 +128,8 @@ classdef Pulseformer % data_out = data_out(st:en); %scaling?! see pulsef module line 696 - scale = max(max([abs(real(data_out)) abs(imag(data_out))])); %find max value from real and imag part - data_out = data_out./scale; +% scale = max(max([abs(real(data_out)) abs(imag(data_out))])); %find max value from real and imag part +% data_out = data_out./scale; data_out = data_out'; diff --git a/Classes/04_DSP/EQ_silas.m b/Classes/04_DSP/EQ_silas.m index 5cf4eb7..d20f438 100644 --- a/Classes/04_DSP/EQ_silas.m +++ b/Classes/04_DSP/EQ_silas.m @@ -12,6 +12,7 @@ classdef EQ_silas < handle d_constellation %constellation points of the reference y_out %equalizer output signal + d_out %decision output % FFE coefficients always named with "e" Ne @@ -107,7 +108,7 @@ classdef EQ_silas < handle end - function [signalclass_out] = process(obj,signalclass_in, reference_signalclass_in) + function [signalclass_out,symbols_out] = process(obj,signalclass_in, reference_signalclass_in) % actual processing of the signal (steps 1. - 3.) % 1 normalize RMS @@ -117,6 +118,7 @@ classdef EQ_silas < handle obj.process_(signalclass_in.signal', reference_signalclass_in.signal'); signalclass_in.signal = obj.y_out'; + %change sampling frequency of outgoing signal signalclass_in.fs = reference_signalclass_in.fs; @@ -125,6 +127,9 @@ classdef EQ_silas < handle lbdesc = ['EQ von Silas ist gelaufen ']; signalclass_in = signalclass_in.logbookentry(lbdesc); + symbols_out = signalclass_in; + symbols_out.signal = obj.d_out; + % write to output signalclass_out = signalclass_in; @@ -213,7 +218,7 @@ classdef EQ_silas < handle %start the dd mode with coefficients from training coeff = [obj.e;obj.b]; - obj.e_dc = ones(obj.eq_updatelatency,1);%.*obj.e_dc; + obj.e_dc = ones(obj.eq_updatelatency,1).*obj.e_dc; dc_block = ones(obj.eq_parallelization_blocklength,1); for ddloop = 1:obj.ddloops @@ -236,8 +241,9 @@ classdef EQ_silas < handle d_feedback = zeros(obj.Cb(1),1); d_vnle = obj.calcVNLENonlinVecs(d_feedback,obj.Ib2,obj.Ib3,obj.Nb,obj.d_norm); d_hat = zeros(obj.x_length,1); - - m_reg = 0; + lvl_err = NaN(obj.x_length,numel(obj.d_constellation)); + lvl_err_mov = NaN(100,numel(obj.d_constellation)); + m_reg = 0; if obj.eq_avg_blocklength > 0 averaging_window = zeros(obj.eq_avg_blocklength,1); @@ -264,19 +270,30 @@ classdef EQ_silas < handle x_d = [x_vnle;-d_vnle]; %Apply filter - %y(m) = (m_reg(end)*dc_cnt + obj.e_dc(end)) + x_d.'* coeff; if obj.mu_dc_dd > 0 y(m) = obj.e_dc(end) + x_d.'* coeff; else y(m) = x_d.'* coeff; end - %Decision + + %Decision 1 + [~,symbol_idx] = min(abs(y(m) - obj.d_constellation)); % decision for closest constellation point + d_hat(k) = obj.d_constellation(symbol_idx); + + %Error between FFE & DFE filtered signal and Decision + obj.error(k) = y(m) - d_hat(k); +% + lvl_err(k,symbol_idx) = obj.error(k); + + lvl_err_mov(:,symbol_idx) = circshift(lvl_err_mov(:,symbol_idx),1); + lvl_err_mov(1,symbol_idx) = obj.error(k); + + %Decision 2 + y(m) = y(m)-mean(lvl_err_mov(:,symbol_idx),'omitnan'); [~,symbol_idx] = min(abs(y(m) - obj.d_constellation)); % decision for closest constellation point d_hat(k) = obj.d_constellation(symbol_idx); - %Error between FFE & DFE filtered signal and Decision - obj.error(k) = y(m) - d_hat(k); %Update FFE and DFE coefficients coeff = coeff - (mu_mat * (obj.error(k) * conj(x_d))); @@ -324,7 +341,34 @@ classdef EQ_silas < handle end end + %% +% +% b = movmean(lvl_err,[500 500],1,"omitnan"); +% +% figure(11) +% for i = 1:4 +% hold on +% stem(lvl_err(:,i)) +% end + + + + %% + obj.y_out = (circshift( y.' ,-(obj.delay))).'; + obj.d_out = d_hat(1:2:end); +% err = obj.error(1:2:end); +% res = NaN(8,length(err)); +% for lvl = 1:8 +% a = find(obj.d_out==obj.d_constellation(lvl)); +% res(lvl,a) = err(a); +% end +% mean(res,2,"omitnan"); +% +% figure(12) +% scatter(1:length(obj.y_out),obj.y_out,1,'.') +% hold on +% scatter(1:length(obj.d_out),obj.d_out,1,'.') end diff --git a/Datatypes/gain_mode.m b/Datatypes/gain_mode.m index 4d63d8c..c998fc5 100644 --- a/Datatypes/gain_mode.m +++ b/Datatypes/gain_mode.m @@ -1,8 +1,8 @@ classdef gain_mode < int32 enumeration - gain (1) - output_power (2) + gain (1) + output_power (2) end end \ No newline at end of file diff --git a/Datatypes/normalization_mode.m b/Datatypes/normalization_mode.m index 768dd8f..acac5eb 100644 --- a/Datatypes/normalization_mode.m +++ b/Datatypes/normalization_mode.m @@ -3,6 +3,7 @@ classdef normalization_mode < int32 enumeration rms (1) oneone (2) + milliwatt (3) end end \ No newline at end of file diff --git a/Datatypes/power_notation.m b/Datatypes/power_notation.m new file mode 100644 index 0000000..bf31d1e --- /dev/null +++ b/Datatypes/power_notation.m @@ -0,0 +1,10 @@ +classdef power_notation < int32 + + enumeration + mW (1) + W (2) + dBW (3) + dBm (4) + end + +end \ No newline at end of file diff --git a/Datatypes/upsampling_mode.m b/Datatypes/upsampling_mode.m new file mode 100644 index 0000000..9e42660 --- /dev/null +++ b/Datatypes/upsampling_mode.m @@ -0,0 +1,8 @@ +classdef upsampling_mode < int32 + + enumeration + resample (1) + samplehold (2) + end + +end \ No newline at end of file diff --git a/Functions/calc_ber.m b/Functions/calc_ber.m index c5cd69a..8322d88 100644 --- a/Functions/calc_ber.m +++ b/Functions/calc_ber.m @@ -50,7 +50,9 @@ end delta_bits = length(reference) - length(data); - skip_end = max(skip_end,delta_bits); +% skip_end = max(skip_end,delta_bits); + skip_end = delta_bits + skip_end; + reference_ = logical(reference(skipstart+1:end-skip_end,:))'; end diff --git a/Functions/plot_eye.m b/Functions/plot_eye.m new file mode 100644 index 0000000..c4fc406 --- /dev/null +++ b/Functions/plot_eye.m @@ -0,0 +1,40 @@ +function plot_eye(signal, fs, fsym) + + + histpoints = 512; %% verticale resolution + histpoints = floor(histpoints/2)*2+1; %% to have the eye digram centered around one point make the vertical resolution uneven + histpoints_horizontal = 256; %% horizontal resolution + hist_data=zeros(histpoints,histpoints_horizontal ); %% initilize eye diagram + x = real(signal); %% make input signal real + + x = resample(x,fsym*histpoints_horizontal/2,fs); %% down sample to original fsym rate + if mod(length(x),2)==1 %% if the signal lenght is not divisible by 2 (symbols displayed in the eye diagram are 2) remove last symbol + x = x(1:end-1); + end + + %x = resample(x,histpoints_horizontal/2,1); %% reshape signal to high resolution eye diagram in this case 2 symbols are upsampled to 512 values + eye_mat = reshape(x(1:end-mod(length(x),histpoints_horizontal)),histpoints_horizontal,floor(length(x)/histpoints_horizontal)); %% reshape signal into 256 rows each row has the histogram(eye data of all symbols) + + %defult value case + maxA = max(eye_mat(:))*1.1; + minA = min(eye_mat(:))*1.1; + + difference= maxA-minA; + + data_ind_y=round((eye_mat-minA)/difference*(histpoints-1)) +1; + + for n=1:size(data_ind_y,1) + [nn,cc] = hist(data_ind_y(n,:),1:histpoints); + hist_data(:,n)=hist_data(:,n)+nn.'; + end + + figure; + image(hist_data); + colormap(cbrewer2("Reds",92)); + + y = string(linspace(minA,maxA,16)); + yticks(linspace(0,histpoints,16)); + yticklabels(y); + + +end \ No newline at end of file diff --git a/projects/MPI_April/auswertung_1.m b/projects/MPI_April/auswertung_1.m index e3ca3a7..1c8275d 100644 --- a/projects/MPI_April/auswertung_1.m +++ b/projects/MPI_April/auswertung_1.m @@ -2,14 +2,14 @@ -vp = wh.parameter.vp.values(1); +vp = wh.parameter.vp.values(2); vb = wh.parameter.vb.values(1); rop = wh.parameter.rop.values; f=figure(1113); -tiledlayout(2,4) -for sir = [20,36] +tiledlayout(2,3) +for sir = [36,26] for lw = wh.parameter.laser_linewidth.values nexttile cols = linspecer(9); @@ -22,8 +22,7 @@ for sir = [20,36] for pn_key = wh.parameter.pn_key.values curber(end+1,:) = wh.getStoValue('ber',sir,lw,pn_key,vp,bias,rop); rop_meas(end+1,:) = wh.getStoValue('rop_save',sir,lw,pn_key,vp,bias,rop); - - curstd(end+1,:,:) = wh.getStoValue('level_std',sir,lw,pn_key,vp,bias,rop); + plot(rop_meas(end,:)-txpow_meas, curber(end,:) ,'LineStyle',':','Color',cols(cnt,:),'LineWidth',0.1,'Marker','o','MarkerEdgeColor',[1 1 1],'MarkerFaceColor',cols(cnt,:),'HandleVisibility','off','MarkerSize',1); hold on end diff --git a/projects/MPI_April/auswertung_2.m b/projects/MPI_April/auswertung_2.m new file mode 100644 index 0000000..01400b1 --- /dev/null +++ b/projects/MPI_April/auswertung_2.m @@ -0,0 +1,42 @@ + +vp = wh.parameter.vp.values(1); +vb = wh.parameter.vb.values(1); +rop = wh.parameter.rop.values; + + +f=figure(113); +tiledlayout(2,4) +for sir = [20,36] + for lw = wh.parameter.laser_linewidth.values + nexttile + cols = linspecer(9); + cnt = 1; + for bias = wh.parameter.vb.values(1:2:end-1) + curber = []; + curstd = []; + rop_meas = []; + txpow_meas = wh.getStoValue('mod_out_pow',sir,lw,1,vp,bias,rop(1)); + + for pn_key = wh.parameter.pn_key.values + rop_meas(end+1,:) = wh.getStoValue('rop_save',sir,lw,pn_key,vp,bias,rop); + curstd(end+1,:,:) = wh.getStoValue('level_std',sir,lw,pn_key,vp,bias,rop); + end + + std_ = squeeze(mean(curstd,1))'; + symbols = ["square",'o','+','x']; + for l =1:size(std_,1) + plot(mean(rop_meas,1)-txpow_meas, std_(l,:) ,'DisplayName',['Vbias: ',num2str(bias),' V; Pam-lvl: ',num2str(l)],'LineStyle','-','Color',cols(cnt,:),'LineWidth',1,'Marker',symbols(l),'MarkerEdgeColor',[1 1 1],'MarkerFaceColor',cols(cnt,:)); + hold on + end + cnt = cnt+1; + end + + title(['SIR: ',num2str(sir),'; Lw: ',num2str(lw*1e-6),' MHz; Vpeakpeak: ',num2str(2*vp)]) + legend('Location','northeast') + xlabel("measured ROP in dBm") + ylabel("Std") + ylim([0.12,0.27]); + drawnow + end +end + diff --git a/projects/MPI_April/auswertung_3.m b/projects/MPI_April/auswertung_3.m new file mode 100644 index 0000000..871b5cc --- /dev/null +++ b/projects/MPI_April/auswertung_3.m @@ -0,0 +1,89 @@ + + +files = ["imdd_simulation\projects\MPI_April\wh_pam4_dienstag.mat","imdd_simulation\projects\MPI_April\wh_pam6_dienstag.mat"]; +linesstyles = ["-","--"]; + +vp = wh.parameter.vp.values; +vb = wh.parameter.vb.values; +sir = wh.parameter.sir.values(end); % 20 22 24 26 28 30 32 34 36 +rop = wh.parameter.rop.values(end); +lw = wh.parameter.laser_linewidth.values(1); +pn_key = wh.parameter.pn_key.values; +bias = wh.parameter.vb.values(1:end-2); +figure() +tiledlayout(2,3) +cols = linspecer(4); + +for d = 1:2 + load(files(d)); + for lw = wh.parameter.laser_linewidth.values + nexttile + + curber = []; + + for v = 1:numel(wh.parameter.vp.values) + vp_ = wh.parameter.vp.values(v); + for b = 1:numel(bias) + + for k = pn_key + curber(k,b) = wh.getStoValue('ber',sir,lw,k,vp_,bias(b),rop); + end + + end + meanber = mean(curber,1); + hold on + + delta = meanber-curber; + yneg = abs(max(delta)); + ypos = abs(min(delta)); + + + + if d == 1 + curvename = ['PAM4, Vpp:',num2str((vp_/2*2)*100),' %']; + [hl,hp] = boundedline(2-bias,meanber,([yneg;ypos]'),'-o','alpha','Color',cols(v,:),'transparency', 0.1,'linewidth',0.7); + hl.MarkerFaceColor = cols(v,:); + hl.MarkerSize = 3; + hl.LineStyle = '-'; + hl.DisplayName = curvename; + set(hp,'HandleVisibility','off'); + + ho = outlinebounds(hl,hp); + set(ho, 'linestyle', ':', 'color', cols(v,:),'Linewidth',0.5); + set(ho,'HandleVisibility','off'); + % errorbar(2-bias,meanber,yneg,ypos,'LineStyle',linesstyles(d),'Marker','o','Color',cols(v,:),'MarkerSize',4,'MarkerFaceColor',cols(v,:),'MarkerEdgeColor',cols(v,:),'DisplayName',curvename); + elseif d == 2 + curvename = ['PAM6, Vpp:',num2str((vp_/2*2)*100),' %']; + [hl,hp] = boundedline(2-bias,meanber,([yneg;ypos]'),'-o','alpha','Color',cols(v,:),'transparency', 0.1,'linewidth',0.7); + hl.MarkerFaceColor = cols(v,:); + hl.MarkerSize = 3; + hl.LineStyle = '--'; + hl.DisplayName = curvename; + set(hp,'HandleVisibility','off'); + + ho = outlinebounds(hl,hp); + set(ho, 'linestyle', ':', 'color', cols(v,:),'Linewidth',0.5); + set(ho,'HandleVisibility','off'); + + % errorbar(2-bias,meanber,yneg,ypos,'LineStyle',linesstyles(d),'Marker','square','Color',cols(v,:),'MarkerSize',4,'MarkerFaceColor',cols(v,:),'MarkerEdgeColor',cols(v,:),'DisplayName',curvename); + end + + % + % plot(2-bias, meanber,'LineStyle',linesstyles(d),'LineWidth',0.5,'DisplayName',curvename,'Color',cols(v,:)); + % scatter(2-bias,curber,4,'o','MarkerFaceColor',cols(v,:),'MarkerEdgeColor',cols(v,:),'HandleVisibility','off'); + hold on + title(['SIR: ',num2str(sir),'; Lw: ',num2str(lw*1e-6),' MHz']) + set(gca,'YScale','log'); + legend('Location','southwest') + xlabel("Bias in V") + ylabel("BER") + + xlim([0.4,1]); + ylim([1e-4,3e-1]); + + end + yline(3.8e-3,'DisplayName','FEC'); + end + +end + diff --git a/projects/MPI_April/auswertung_4.m b/projects/MPI_April/auswertung_4.m new file mode 100644 index 0000000..56d3e04 --- /dev/null +++ b/projects/MPI_April/auswertung_4.m @@ -0,0 +1,52 @@ + + +files = ["imdd_simulation\projects\MPI_April\wh_pam4_dienstag.mat","imdd_simulation\projects\MPI_April\wh_pam6_dienstag.mat"]; +linesstyles = ["-","--"]; + +vp = wh.parameter.vp.values; +vb = wh.parameter.vb.values; +sir = wh.parameter.sir.values(4); +rop = wh.parameter.rop.values(end); +lw = wh.parameter.laser_linewidth.values(1); +pn_key = wh.parameter.pn_key.values; +bias = wh.parameter.vb.values(1:end-3); + +figure() +tiledlayout(2,3) +cols = flip(cbrewer2("RdYlBu",32)); + +for d = 1:2 + load(files(d)); + + for lw = wh.parameter.laser_linewidth.values + nexttile + + curber = []; + meanber = []; + for v = 1:numel(wh.parameter.vp.values) + vp_ = wh.parameter.vp.values(v); + for b = 1:numel(bias) + for k = pn_key + curber(k,b) = wh.getStoValue('ber',sir,lw,k,vp_,bias(b),rop); + end + end + meanber(v,:) = mean(curber,1); + end + + levels = [1e-4, reshape([1e-4; 1e-3; 1e-2]*[2:2:10],1,[])]; + contourf(2-bias,wh.parameter.vp.values*100, meanber,levels,'DisplayName',"BER","EdgeAlpha",0.5,'LineStyle','--'); + colormap(cols) + set(gca,'ColorScale','log'); + xlabel('bias ') + ylabel('vpp in %') + clim([1e-4 3e-1]); + + if d == 1 + title(['PAM4 - SIR: ',num2str(sir),'; Lw: ',num2str(lw*1e-6),' MHz']) + elseif d == 2 + title(['PAM6 - SIR: ',num2str(sir),'; Lw: ',num2str(lw*1e-6),' MHz']) + end + + end +end + diff --git a/projects/MPI_April/auswertung_5.m b/projects/MPI_April/auswertung_5.m new file mode 100644 index 0000000..ed49c2f --- /dev/null +++ b/projects/MPI_April/auswertung_5.m @@ -0,0 +1,49 @@ + + +files = ["imdd_simulation\projects\MPI_April\wh_pam4_dienstag.mat","imdd_simulation\projects\MPI_April\wh_pam6_dienstag.mat"]; +linesstyles = ["-","--"]; + +vp = wh.parameter.vp.values; +vb = wh.parameter.vb.values; +sir = wh.parameter.sir.values(end); +rop = wh.parameter.rop.values(end); +lw = wh.parameter.laser_linewidth.values(1); +pn_key = wh.parameter.pn_key.values; +bias = wh.parameter.vb.values(1:end-3); + +figure() +tiledlayout(2,3) +cols = flip(cbrewer2("RdYlBu",32)); + +for d = 1:2 + load(files(d)); + nexttile + for lw = wh.parameter.laser_linewidth.values + + curcspr = []; + meancspr = []; + for v = 1:numel(wh.parameter.vp.values) + vp_ = wh.parameter.vp.values(v); + + for b = 1:numel(bias) + for k = pn_key + curcspr(k,b) = wh.getStoValue('cspr',sir,lw,k,vp_,bias(b),rop); + end + end + + meancspr(v,:) = mean(curcspr,1); + hold on + plot(2-bias,meancspr(v,:),'DisplayName',['vpp: ',num2str(vp_)]); + + end + + + if d == 1 + title(['PAM4 - SIR: ',num2str(sir),'; Lw: ',num2str(lw*1e-6),' MHz']) + elseif d == 2 + title(['PAM6 - SIR: ',num2str(sir),'; Lw: ',num2str(lw*1e-6),' MHz']) + end + + end +end + diff --git a/projects/MPI_April/mpi_simulation_cspr.m b/projects/MPI_April/mpi_simulation_cspr.m index e83439f..482055e 100644 --- a/projects/MPI_April/mpi_simulation_cspr.m +++ b/projects/MPI_April/mpi_simulation_cspr.m @@ -1,214 +1,236 @@ %% Settings + clear +for M=[4,6,8] -filename = '112G_2'; -load_sequence = 0; -M = 6; -datarate = 448e9; + filename = '112G_2'; + load_sequence = 0; -kover = 4; -fsym = round(datarate*1e-9 / log2(M))*1e9; -fdac = 256e9; -fadc = 256e9; + datarate = 448e9; -lowpass_cutoff = fsym/2 * 1.1; -awg_bw = lowpass_cutoff; -mod_bw = lowpass_cutoff; -phd_bw = lowpass_cutoff; -scp_bw = lowpass_cutoff; + kover = 8; + fsym = round(datarate*1e-9 / log2(M))*1e9; + fdac = 256e9; + fadc = 256e9; -LP_awg = Filter('filtdegree',4,"f_cutoff",90e9,"fs",fdac*kover,"filterType",filtertypes.butterworth); -LP_modulator= Filter('filtdegree',2,"f_cutoff",70e9,"fs",fdac*kover,"filterType",filtertypes.butterworth); -LP_opt = Filter('filtdegree',3,"f_cutoff",fsym/log2(M).*1.5,"fs",fdac*kover,"filterType",filtertypes.gaussian); -LP_phd = Filter('filtdegree',2,"f_cutoff",70e9,"fs",fdac*kover,"filterType",filtertypes.butterworth); -LP_scpe = Filter('filtdegree',4,"f_cutoff",110e9,"fs",fadc,"filterType",filtertypes.butterworth); + lowpass_cutoff = fsym/2 * 1.1; + awg_bw = lowpass_cutoff; + mod_bw = lowpass_cutoff; + phd_bw = lowpass_cutoff; + scp_bw = lowpass_cutoff; + + LP_awg = Filter('filtdegree',4,"f_cutoff",90e9,"fs",fdac*kover,"filterType",filtertypes.butterworth); + LP_modulator= Filter('filtdegree',2,"f_cutoff",70e9,"fs",fdac*kover,"filterType",filtertypes.butterworth); + LP_opt = Filter('filtdegree',3,"f_cutoff",fsym/log2(M).*1.5,"fs",fdac*kover,"filterType",filtertypes.gaussian); + LP_phd = Filter('filtdegree',2,"f_cutoff",70e9,"fs",fdac*kover,"filterType",filtertypes.butterworth); + LP_scpe = Filter('filtdegree',4,"f_cutoff",110e9,"fs",fadc,"filterType",filtertypes.butterworth); + + + % 1) PRBS Generation + O = 18; %order of prbs + N = 2^(O-1); %length of prbs + [~,seed] = prbs(O,1); %initialize first seed of prbs + bitpattern=[]; + for i = 1:log2(M) + [bitpattern(:,i),seed] = prbs(O,N,seed); + end + if M == 6 + bitpattern = reshape(bitpattern,[],1); + bitpattern = bitpattern(1:end-mod(length(bitpattern),5)); + end + + % 2 ) Build Inf. signal class + bits = Informationsignal(bitpattern); + + % 3) Digi modulation -> PAM-M signal + digimod_out = PAMmapper(M,0).map(bits); + digimod_out.fs = fsym; + + + % linearGain = 1; + % limit = 1; + % SaturatingAmplifier = serdes.SaturatingAmplifier('Mode',1,... + % 'Limit',limit,'LinearGain',linearGain); + % X.signal = SaturatingAmplifier(X.signal); + % X.signal = min(max(X.signal,-0.8),0.8); + % X = X.normalize("mode","oneone"); + + sir = [20:2:36]; %decibel = attenuation of interference path + laser_linewidth = [1e5 1e6 10e6]; + pn_key = [1:10]; + vp = [0.25,0.5,0.75,1]; + vb = [1:0.1:1.8]; + + rop = -5:5; + sir = 28; + laser_linewidth = 10e6; + pn_key = 9; + vp = 0.5; + vb = 1;%[1:0.1:1.8]; + mpi_path=50; + + cnt = 1; + + for s = 1:length(sir) + for l = 1:length(laser_linewidth) + for pnk = 1:length(pn_key) + for n = 1:length(vp) + for m = 1:length(vb) + + %digimod_out = digimod_out.normalize("mode","oneone"); + + % cnt = cnt+1; + X = Pulseformer("fsym",fsym,"fdac",fdac,"pulse","rrc","pulselength",16,"rrcalpha",0.5).process(digimod_out); + + % 5) AWG (lowpass, quantization, sample and hold) + X = AWG("fdac",fdac,"dac_min",-1,"dac_max",1,"H_lpf",LP_awg,"kover",kover,"bit_resolution",5,"lpf_active",1).process(X); + + % 6) Lowpass behavior before laser + X = LP_modulator.process(X); + + % 7) Normalize signal + X = X.normalize("mode","oneone"); + + % 1) Laser; Modulation -> OPTICAL DOMAIN + u_pi = 2; + vbias = -vb(m); + extmodlaser = EML("mode",eml_mode.im_cosinus,"power",3,"fsimu",X.fs,"lambda",1290,"bias",vbias,"u_pi",u_pi,"linewidth",laser_linewidth(l),"randomkey",pn_key(pnk)); + E = X.*vp(n); + + [Opt,extmodlaser] = extmodlaser.process(E); + + % figure(m) + % hold on + % scatter(E.signal(1:100000),(abs(Opt.signal(1:100000)).^2)*1e3,0.1,'.','DisplayName','Modulator TF') + % xlabel('Input in V') + % ylabel('abs(Output) in mW') + + % ER = 10*log10(max(abs(Opt.signal).^2)/min(abs(Opt.signal).^2)); + + Opt = LP_opt.process(Opt); + + cspr(s,l,pnk,n,m) = Opt.cspr; + mod_out_pow(s,l,pnk,n,m) = Opt.power; + + % 2) ping pong fiber propagation + Interference_sig = Fiber("fsimu",Opt.fs,"fiber_length",mpi_path*2/1000,"alpha",0,"D",0,"lambda0",1310,"gamma",0).process(Opt); + + Interference_sig = Amplifier("amp_mode","ideal_no_noise","gain_mode","gain","amplification_db",-sir(s)).process(Interference_sig); + + % In the meantime: delay the main signal + [Main_sig,dly] = Opt.delay("delay_meter",mpi_path*2); + + % Add + Combined_sig = Main_sig + Interference_sig; + + % Cut (due to the delays there is a jump in the signals) + if dly == 0;dly = 1;end + Combined_sig.signal = Combined_sig.signal(ceil(dly):end); + + % Fiber + Combined_sig = Fiber("fsimu",Combined_sig.fs,"fiber_length",2,"alpha",0.3,"D",0,"lambda0",1310,"gamma",0,"Dslope",0.08).process(Combined_sig); + + parfor i = 1:length(rop) + + % Set ROP + Rx_sig = Amplifier("amp_mode","ideal_no_noise","gain_mode","gain","amplification_db",rop(i)).process(Combined_sig); + rop_save(s,l,pnk,n,m,i) = Rx_sig.power; + + % Square Law + Rx_sig = Photodiode("fsimu",Rx_sig.fs,"dark_current",2e-08,"responsivity",1,"temperature",20).process(Rx_sig); + + %Lowpass PhDiode + Rx_sig = LP_phd.process(Rx_sig); + + % Scope + Scpe_sig = Scope("fsimu",Rx_sig.fs,"fadc",fadc,... + "delay",0,"fixed_delay",0,"lpf_bw",scp_bw,"filtertype",filtertypes.butterworth,... + "samplingdelay",0,"rand_samplingdelay",0,"freq_offset",0,"samp_jitter",0,... + "adcresolution",16,"quantbuffer",0.1,'block_dc',1,'lpf_active',1,'H_lpf',LP_scpe).process(Rx_sig); + + % Sample to 2x fsym + Scpe_sig = Scpe_sig.resample("fs_in",fadc,"fs_out",2*fsym); + + % Sync Rx signal with reference + [Scpe_sig,D,cuts] = Scpe_sig.tsynch("reference",digimod_out,"fs_ref",fsym); + + % % % simple EQ (optimum mudc: 0.05 -> 0.005) + % EQ_sig = EQ_silas_plain("Ne",[20,8,8],"Nb",[2,0,0],"trainlength",4096,"mu_dc_dd",0.005,"mu_dc_train",0.05,... + % "mu_ffe_train",0.005,"mu_combined_dd",[0.0004 0.0006 0.0003 0.005],"ddloops",3,'trainloops',3,'sps',2).process(Rx_sig,digimod_out); + +% EQ_sig = EQ("K",2,"plottrain",0,"plotfinal",0,... +% "training_length",4096,"training_loops",3,... +% "Ne",[50,8,8],"Nb",[2,0,0],... +% "DCmu",0.005,"DDmu",[0.0004 0.0006 0.0003 0.005],"DFEmu",0.005,"FFEmu",0.00,... +% "dd_loops",3,"epsilon",[10 100 1000 ],"M",2,... +% "thres",[0.005 0.004 0.0005 ],"l1act",0,"delay",0,"rho",0.0005,"ideal_dfe",0,"DB_aim",0).process(Scpe_sig,digimod_out); + + [EQ_sig,EQ_sym] = EQ_silas("Ne",[50,8,8],"Nb",[2,0,0],"trainlength",4096,... + "sps",2,... + "mu_dc_dd",0.00,... + "mu_dc_train",0.0,... + "mu_ffe_train",0.00,... + "mu_dfe_train",0.005,... + "mu_ffe_dd",[0.0004 0.0006 0.0003],... + "mu_dfe_dd",0.005,... + "ddloops",3,... + "trainloops",3,... + "eq_parallelization_blocklength",1, ... + "eq_updatelatency",1,... + "eq_avg_blocklength",0).process(Scpe_sig,digimod_out); + + % Demap +% Rx_Bits = PAMmapper(M,0).demap(EQ_sig); + Rx_Bits = PAMmapper(M,0).demap(EQ_sym); + % + % Rx_symboldecision = PAMmapper(M,0).decide_pamlevel(EQ_sig,"symbol_levels",unique(digimod_out.signal)); + % + % levels = PAMmapper(M,0).separate_pamlevels(EQ_sig); + % %levels = PAMmapper(M,0).separate_pamlevels(Rx_sig.resample("fs_in",Rx_sig.fs,"fs_out",fsym)); + % + % level_avg(s,l,pnk,n,m,i,:) = mean(levels,'omitnan'); + % level_std(s,l,pnk,n,m,i,:) = std(levels,'omitnan'); + + % BER + [~,errors_bm,BER(s,l,pnk,n,m,i),errors] = calc_ber(Rx_Bits.signal,bitpattern,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + formatted_ber = sprintf('%.1e', BER(s,l,pnk,n,m,i)); + disp(['SIR: ',num2str(sir(s)),'; Lw:',num2str(laser_linewidth(l)),'; Key:',num2str(pn_key(pnk)),'; Vpeak: ',num2str(vp(n)),'; Vbias',num2str(vbias),'; BER: ',formatted_ber,'; run: ',num2str(cnt),' / 12961']); + + % plot_analysis_window; + % drawnow; + + end + end + end + save('pam4_level_comp'); + disp('saved_run2'); + end + end + end + + + % BER plot + figure(340) + cols = linspecer(7); + for m = 1:size(BER,1) + hold on + plot(rop,BER(m,:),'DisplayName',['Bias: ',num2str(vb(m)), ' V; PAM', num2str(M)],'LineStyle','--','Color',cols(M/2,:),'LineWidth',1,'Marker','square','MarkerEdgeColor',[1 1 1],'MarkerFaceColor',cols(M/2,:)); + end + set(gca,'YScale','log'); + legend + xlabel("ROP in dBm") + yline(3.8e-3,'DisplayName','FEC'); + +end figure(2) LP_awg.showHere; -% LP_laser.showHere; -% LP_opt.showHere; -% LP_phd.showHere; -% LP_scpe.showHere; - -% 1) PRBS Generation -O = 18; %order of prbs -N = 2^(O-1); %length of prbs -[~,seed] = prbs(O,1); %initialize first seed of prbs -bitpattern=[]; -for i = 1:log2(M) - [bitpattern(:,i),seed] = prbs(O,N,seed); -end -if M == 6 - bitpattern = reshape(bitpattern,[],1); - bitpattern = bitpattern(1:end-mod(length(bitpattern),5)); -end - -% 2 ) Build Inf. signal class -bits = Informationsignal(bitpattern); - -% 3) Digi modulation -> PAM-M signal -digimod_out = PAMmapper(M,0).map(bits); -digimod_out.fs = fsym; - -X = Pulseformer("fsym",fsym,"fdac",fdac,"pulse","rrc","pulselength",16,"rrcalpha",0.1).process(digimod_out); - -% 5) AWG (lowpass, quantization, sample and hold) -X = AWG("fdac",fdac,"dac_min",-1,"dac_max",1,"H_lpf",LP_awg,"kover",kover,"bit_resolution",5,"lpf_active",1).process(X); - -% 6) Lowpass behavior before laser -X = LP_modulator.process(X); - -% 7) Normalize signal -X = X.normalize("mode","oneone"); -% linearGain = 1; -% limit = 1; -% SaturatingAmplifier = serdes.SaturatingAmplifier('Mode',1,... -% 'Limit',limit,'LinearGain',linearGain); -% X.signal = SaturatingAmplifier(X.signal); -% X.signal = min(max(X.signal,-0.8),0.8); -% X = X.normalize("mode","oneone"); - -sir = [20:2:36]; %decibel = attenuation of interference path - laser_linewidth = [1e5 1e6 10e6]; - pn_key = [1:10]; - vp = [0.25,0.5,0.75,1]; - vb = [1:0.1:1.8]; - - rop = -5:0; -% sir = 35; -% laser_linewidth = 1e6; -% pn_key = 9; -% vp = 0.5; -% vb = [1:0.1:1.8]; -mpi_path=50; - -cnt = 1; -for s = 1:length(sir) - for l = 1:length(laser_linewidth) - for pnk = 1:length(pn_key) - for n = 1:length(vp) - parfor m = 1:length(vb) - - % cnt = cnt+1; - - % 1) Laser; Modulation -> OPTICAL DOMAIN - u_pi = 2; - vbias = -vb(m); - extmodlaser = EML("mode",eml_mode.im_cosinus,"power",3,"fsimu",X.fs,"lambda",1290,"bias",vbias,"u_pi",u_pi,"linewidth",laser_linewidth(l),"randomkey",pn_key(pnk)); - E = X.*vp(n); - - [Opt,extmodlaser] = extmodlaser.process(E); - - % figure(m) - % hold on - % scatter(E.signal(1:100000),(abs(Opt.signal(1:100000)).^2)*1e3,0.1,'.','DisplayName','Modulator TF') - % xlabel('Input in V') - % ylabel('abs(Output) in mW') - - - % ER = 10*log10(max(abs(Opt.signal).^2)/min(abs(Opt.signal).^2)); - - Opt = LP_opt.process(Opt); - - cspr(s,l,pnk,n,m) = Opt.cspr; - mod_out_pow(s,l,pnk,n,m) = Opt.power; - - % 2) ping pong fiber propagation - Interference_sig = Fiber("fsimu",Opt.fs,"fiber_length",mpi_path*2/1000,"alpha",0,"D",0,"lambda0",1310,"gamma",0).process(Opt); - - Interference_sig = Amplifier("amp_mode","ideal_no_noise","gain_mode","gain","amplification_db",-sir(s)).process(Interference_sig); - - % In the meantime: delay the main signal - [Main_sig,dly] = Opt.delay("delay_meter",mpi_path*2); - - % Add - Combined_sig = Main_sig + Interference_sig; - - % Cut (due to the delays there is a jump in the signals) - if dly == 0;dly = 1;end - Combined_sig.signal = Combined_sig.signal(ceil(dly):end); - - % Fiber - Combined_sig = Fiber("fsimu",Combined_sig.fs,"fiber_length",2,"alpha",0.3,"D",0,"lambda0",1310,"gamma",0,"Dslope",0.08).process(Combined_sig); - - for i = 1:length(rop) - - % Set ROP - Rx_sig = Amplifier("amp_mode","ideal_no_noise","gain_mode","gain","amplification_db",rop(i)).process(Combined_sig); - rop_save(s,l,pnk,n,m,i) = Rx_sig.power; - - % Square Law - Rx_sig = Photodiode("fsimu",Rx_sig.fs,"dark_current",2e-08,"responsivity",1,"temperature",20).process(Rx_sig); - - %Lowpass PhDiode - Rx_sig = LP_phd.process(Rx_sig); - - % Scope - Scpe_sig = Scope("fsimu",Rx_sig.fs,"fadc",fadc,... - "delay",0,"fixed_delay",0,"lpf_bw",scp_bw,"filtertype",filtertypes.butterworth,... - "samplingdelay",0,"rand_samplingdelay",0,"freq_offset",0,"samp_jitter",0,... - "adcresolution",16,"quantbuffer",0.1,'block_dc',1,'lpf_active',1,'H_lpf',LP_scpe).process(Rx_sig); - - % Sample to 2x fsym - Scpe_sig = Scpe_sig.resample("fs_in",fadc,"fs_out",2*fsym); - - % Sync Rx signal with reference - [Scpe_sig,D,cuts] = Scpe_sig.tsynch("reference",digimod_out,"fs_ref",fsym); - - % % % simple EQ (optimum mudc: 0.05 -> 0.005) - % EQ_sig = EQ_silas_plain("Ne",[20,8,8],"Nb",[2,0,0],"trainlength",4096,"mu_dc_dd",0.005,"mu_dc_train",0.05,... - % "mu_ffe_train",0.005,"mu_combined_dd",[0.0004 0.0006 0.0003 0.005],"ddloops",3,'trainloops',3,'sps',2).process(Rx_sig,digimod_out); - % -% EQ_sig = EQ("K",2,"plottrain",0,"plotfinal",0,... -% "training_length",4096,"training_loops",3,... -% "Ne",[50,8,8],"Nb",[2,0,0],... -% "DCmu",0.00,"DDmu",[0.0004 0.0006 0.0003 0.005],"DFEmu",0.005,"FFEmu",0.00,... -% "dd_loops",3,"epsilon",[10 100 1000 ],"M",2,... -% "thres",[0.005 0.004 0.0005 ],"l1act",0,"delay",0,"rho",0.0005,"ideal_dfe",0,"DB_aim",0).process(Scpe_sig,digimod_out); - - EQ_sig = EQ_silas("Ne",[50,8,8],"Nb",[2,0,0],"trainlength",4096,... - "sps",2,... - "mu_dc_dd",0.00,... - "mu_dc_train",0.00,... - "mu_ffe_train",0.00,... - "mu_dfe_train",0.005,... - "mu_ffe_dd",[0.0004 0.0006 0.0003],... - "mu_dfe_dd",0.005,... - "ddloops",3,... - "trainloops",3,... - "eq_parallelization_blocklength",1, ... - "eq_updatelatency",1,... - "eq_avg_blocklength",0).process(Scpe_sig,digimod_out); - - % Demap - Rx_Bits = PAMmapper(M,0).demap(EQ_sig); -% -% Rx_symboldecision = PAMmapper(M,0).decide_pamlevel(EQ_sig,"symbol_levels",unique(digimod_out.signal)); -% -% levels = PAMmapper(M,0).separate_pamlevels(EQ_sig); -% %levels = PAMmapper(M,0).separate_pamlevels(Rx_sig.resample("fs_in",Rx_sig.fs,"fs_out",fsym)); -% -% level_avg(s,l,pnk,n,m,i,:) = mean(levels,'omitnan'); -% level_std(s,l,pnk,n,m,i,:) = std(levels,'omitnan'); - - % BER - [~,errors_bm,BER(s,l,pnk,n,m,i),errors] = calc_ber(Rx_Bits.signal,bitpattern,"skip_front",0,"skip_end",0,"returnErrorLocation",1); - - formatted_ber = sprintf('%.1e', BER(s,l,pnk,n,m,i)); - disp(['SIR: ',num2str(sir(s)),'; Lw:',num2str(laser_linewidth(l)),'; Key:',num2str(pn_key(pnk)),'; Vpeak: ',num2str(vp(n)),'; Vbias',num2str(vbias),'; BER: ',formatted_ber,'; run: ',num2str(cnt),' / 12961']); - -% plot_analysis_window; -% drawnow; - - end - end - end - save('pam6_level_comp'); - disp('saved_run2'); - end - end -end +LP_modulator.showHere; +LP_opt.showHere; +LP_phd.showHere; +LP_scpe.showHere; % figure('Name','spectrum') % tiledlayout(4,1) @@ -224,21 +246,9 @@ end % save(['C:\Users\Silas\Nextcloud4\Dokumente\02_Ablage_Office\MPI\Investigation_April_2024\','PAM_',num2str(M),'_mpi_',num2str(mpi_path),'_lw_',num2str(laser_linewidth)],"BER"); -% -% % BER plot -% figure(340) -% cols = linspecer(7); -% for m = 1:size(BER,1) -% hold on -% plot(rop,BER(m,:),'DisplayName',['Bias: ',num2str(vb(m)), ' V'],'LineStyle','-','Color',cols(m,:),'LineWidth',1,'Marker','o','MarkerEdgeColor',[1 1 1],'MarkerFaceColor',cols(m,:)); -% -% end -% set(gca,'YScale','log'); -% legend -% xlabel("ROP in dBm") -% yline(3.8e-3,'DisplayName','FEC'); -% -% + +% +% % figure(21) % hold on % plot(rop,mean(BER),'DisplayName',['Modulation: ',num2str(2*vp/extmodlaser.u_pi*100), ' $\%$'],'LineStyle','-','Color',cols(2,:),'LineWidth',1); @@ -246,14 +256,14 @@ end % legend % ylabel("ROP in dBm") % yline(3.8e-3,'DisplayName','FEC'); -% -% -% -% -% -% -% -% +% +% +% +% +% +% +% +% % %check Rx and TX symbols % figure(101) % scatter(1:100,Rx_symboldecision.signal(1:100),10,'o'); diff --git a/projects/MPI_April/plot_analysis_window.m b/projects/MPI_April/plot_analysis_window.m index 626938b..4f1b137 100644 --- a/projects/MPI_April/plot_analysis_window.m +++ b/projects/MPI_April/plot_analysis_window.m @@ -16,9 +16,9 @@ ylim([-100,0]); %Phase Investigation subplot(3,2,3) hold on -phase_int = extmodlaser.phase(ceil(n):end); -phase_main = delayseq(extmodlaser.phase,ceil(n)); -phase_main = phase_main(ceil(n):end); +phase_int = extmodlaser.phase(ceil(dly):end); +phase_main = delayseq(extmodlaser.phase,ceil(dly)); +phase_main = phase_main(ceil(dly):end); phase_diff = phase_int - phase_main; t = (1:length(phase_int))' ./ extmodlaser.fsimu ; plot(t*1e6,phase_int,'Color',colpairs(1,:),'DisplayName','Interferer Phase'); @@ -30,7 +30,7 @@ subplot(3,2,4) %Received Signal after Phdiode yyaxis left t = (1:Rx_sig.length)' ./ Rx_sig.fs; -scatter(t*1e6,Rx_sig.normalize("mode","oneone").signal.*max(unique(digimod_out.signal)),1,'.','MarkerEdgeColor',cols(5,:)); +scatter(t*1e6,Rx_sig.normalize("mode","oneone").signal.*max(unique(digimod_out.signal)),1,'.','MarkerEdgeColor',cols(6,:)); hold on errors_t = errors./Rx_Bits.fs; errors_sym = Rx_symboldecision.signal(errors); @@ -54,7 +54,7 @@ for h = 1:size(levels,2) hold on end %yline(PAMmapper(M,0).thresholds); -title(['BER: ',sprintf('%.1e', BER(m,i))]); +title(['BER: ',sprintf('%.1e', BER(s,l,pnk,n,m,i))]); %Histogram of EQzed Signal subplot(3,2,6) diff --git a/projects/MPI_April/save_to_warehouse.m b/projects/MPI_April/save_to_warehouse.m index f79e282..42bf10d 100644 --- a/projects/MPI_April/save_to_warehouse.m +++ b/projects/MPI_April/save_to_warehouse.m @@ -6,13 +6,20 @@ vp = [0.25]; vb = [1:0.1:1.8]; rop = -9:3; +sir = [20:2:36]; %decibel = attenuation of interference path +laser_linewidth = [1e5 1e6 10e6]; +pn_key = [1:10]; +vp = [0.25,0.5,0.75,1]; +vb = [1:0.1:1.8]; +rop = -5:0; + params = struct; params.sir = [20:2:36]; %decibel = attenuation of interference path -params.laser_linewidth = [1e5 1e6 3e6 10e6]; +params.laser_linewidth = [1e5 1e6 10e6]; params.pn_key = [1:10]; params.vp = [0.25,0.5,0.75,1]; params.vb = [1:0.1:1.8]; -params.rop = -9:3; +params.rop = -5:0; wh = DataStorage(params); wh.addStorage("ber"); wh.addStorage("level_avg"); @@ -21,18 +28,19 @@ wh.addStorage("rop_save"); wh.addStorage("cspr"); wh.addStorage("mod_out_pow"); + cnt = 1; for s = 1:length(sir) for l = 1:length(laser_linewidth) for pnk = 1:length(pn_key) - for n = 1 + for n = 1:length(vp) for m = 1:length(vb) for i = 1:length(rop) cnt = cnt +1; wh.addValueToStorage(BER(s,l,pnk,n,m,i),'ber',sir(s),laser_linewidth(l),pn_key(pnk),vp(n),vb(m),rop(i)); - wh.addValueToStorage(level_avg(s,l,pnk,n,m,i,:),'level_avg',sir(s),laser_linewidth(l),pn_key(pnk),vp(n),vb(m),rop(i)); - wh.addValueToStorage(level_std(s,l,pnk,n,m,i,:),'level_std',sir(s),laser_linewidth(l),pn_key(pnk),vp(n),vb(m),rop(i)); +% wh.addValueToStorage(level_avg(s,l,pnk,n,m,i,:),'level_avg',sir(s),laser_linewidth(l),pn_key(pnk),vp(n),vb(m),rop(i)); +% wh.addValueToStorage(level_std(s,l,pnk,n,m,i,:),'level_std',sir(s),laser_linewidth(l),pn_key(pnk),vp(n),vb(m),rop(i)); wh.addValueToStorage(rop_save(s,l,pnk,n,m,i),'rop_save',sir(s),laser_linewidth(l),pn_key(pnk),vp(n),vb(m),rop(i)); wh.addValueToStorage(cspr(s,l,pnk,n,m),'cspr',sir(s),laser_linewidth(l),pn_key(pnk),vp(n),vb(m),rop(i)); wh.addValueToStorage(mod_out_pow(s,l,pnk,n,m),'mod_out_pow',sir(s),laser_linewidth(l),pn_key(pnk),vp(n),vb(m),rop(i)); diff --git a/projects/MPI_April/wh_pam4_dienstag.mat b/projects/MPI_April/wh_pam4_dienstag.mat new file mode 100644 index 0000000..ce6b642 Binary files /dev/null and b/projects/MPI_April/wh_pam4_dienstag.mat differ diff --git a/projects/MPI_April/wh_pam6_dienstag.mat b/projects/MPI_April/wh_pam6_dienstag.mat new file mode 100644 index 0000000..e463b5b Binary files /dev/null and b/projects/MPI_April/wh_pam6_dienstag.mat differ diff --git a/test/modulator_test.m b/test/modulator_test.m new file mode 100644 index 0000000..8ae7ffe --- /dev/null +++ b/test/modulator_test.m @@ -0,0 +1,105 @@ + +datarate = 112e9; +M = 6; +laser_linewidth = 0; +kover = 32; +fsym = round(datarate*1e-9 / log2(M))*1e9; +fdac = 256e9; + +% 1) PRBS Generation +O = 18; %order of prbs +N = 2^(O-1); %length of prbs +[~,seed] = prbs(O,1); %initialize first seed of prbs +bitpattern=[]; + +for i = 1:log2(M) + [bitpattern(:,i),seed] = prbs(O,N,seed); +end + +if M == 6 + bitpattern = reshape(bitpattern,[],1); + bitpattern = bitpattern(1:end-mod(length(bitpattern),5)); +end + +% 2 ) Build Inf. signal class +bits = Informationsignal(bitpattern); + + +alphas = flip([0.01,0.05,0.1,0.2,0.3,0.4]); +LP_awg = Filter('filtdegree',5,"f_cutoff",70e9,"fs",fdac*kover,"filterType",filtertypes.butterworth); +AWG_=AWG("fdac",fdac,"dac_min",-1,"dac_max",1,"H_lpf",LP_awg,"kover",kover,"bit_resolution",16,"lpf_active",1,"normalize2dac",1,"upsampling_method","samplehold"); + + +% 3) Digi modulation -> PAM-M signal +digimod_out = PAMmapper(M,0).map(bits); +digimod_out.fs = fsym; + +X = Pulseformer("fsym",fsym,"fdac",fdac,"pulse","rrc","pulselength",16,"rrcalpha",alphas(a)).process(digimod_out); + +% plot_eye(X.signal,X.fs,fsym); +% 5) AWG (lowpass, quantization, sample and hold) +X = AWG_.process(X); + +f=figure(12); +f.Name = 'spectrum'; +spectrum_plot(X.signal,X.fs,'spectrum'); + +figure(4) +LP_awg.showHere(); + +disp(['El. power: ',num2str(X.power),' dBm (into 50 Ohm)']); +disp(['El. RMS voltage: ',num2str(sqrt(mean(X.signal.^2))),' V']); +disp(['max voltage: ',num2str(max(X.signal)),' V']); + +rms_ = rms(X.signal); +max_ = max(X.signal); +min_ = min(X.signal); +figure(1); +hold on +plot(X.signal,'LineWidth',0.1); +yline([max_, min_],'LineWidth',2,'LineStyle','--'); +yline([rms_, -rms_],'LineWidth',2,'LineStyle',':'); +ylim([-3 3]); +title(['AWG output: ',num2str(X.power), 'dBm']); + +f=figure(12); +f.Name = 'spectrum'; +spectrum_plot(X.signal,X.fs,'spectrum'); + +% 6) Lowpass behavior before laser +LP_modulator= Filter('filtdegree',4,"f_cutoff",70e9,"fs",fdac*kover,"filterType",filtertypes.butterworth); +X = LP_modulator.process(X); + + +% 7) Normalize signal to 0dB rms +figure(1);subplot(1,2,1);plot(X.signal,'LineWidth',0.1);ylim([-3 3]);title(['AWG output: ',num2str(X.power), 'dBm into 50 Ohm']); + +% 1) Laser; Modulation -> OPTICAL DOMAIN +u_pi = 2; +vbias = -1; +extmodlaser = EML("mode",eml_mode.im_cosinus,"power",0,"fsimu",X.fs,"lambda",1290,"bias",vbias,"u_pi",u_pi,"linewidth",laser_linewidth,"randomkey",5); + +E = X.normalize("mode","oneone"); +figure(1);subplot(1,2,2);plot(E.signal,'LineWidth',0.1);ylim([-3 3]);title(['scaled to modulator; ',num2str(E.power), 'dBm']); +disp(['El. power: ',num2str(E.power),' dBm (into 50 Ohm)']); + +[Opt,extmodlaser] = extmodlaser.process(E); + +f=figure(12); +f.Name = 'spectrum'; +spectrum_plot(X.signal,X.fs,'spectrum'); + +plot_eye(E.signal,E.fs,fsym); + +figure(111) +hold on +scatter(E.signal(1:100000),(abs(Opt.signal(1:100000)).^2)*1e3,0.1,'.','DisplayName','Modulator TF') +xlabel('Input in V') +ylabel('Output in mW') + + + + + + + diff --git a/test/test_modulation.m b/test/test_modulation.m new file mode 100644 index 0000000..e080981 --- /dev/null +++ b/test/test_modulation.m @@ -0,0 +1,36 @@ + + + + +O = 10; %order of prbs +N = 2^(O-1); %length of prbs +[~,seed] = prbs(O,1); %initialize first seed of prbs +bitpattern=[]; +M=6; + +for i = 1:log2(M) + [bitpattern(:,i),seed] = prbs(O,N,seed); +end + +if M == 6 + bitpattern = reshape(bitpattern,[],1); +end + +% [bitpattern,seed] = prbs(O,N,seed); +% % +% bitpattern = bitpattern'; + +% 2 ) Build Inf. signal class +bits = Informationsignal(bitpattern); + + +% 3) Digi modulation -> PAM-M signal +digimod_out = PAMmapper(M,0).map(bits); + +Rx_Bits = PAMmapper(M,0).demap(digimod_out); + +[~,errors_bm,BER,errors] = calc_ber(Rx_Bits.signal,bitpattern,"skip_front",0,"skip_end",0,"returnErrorLocation",1); + +formatted_ber = sprintf('%.1e', BER); +disp(['BER: ',formatted_ber]); + \ No newline at end of file