diff --git a/Time-Series-Analyzer/computeAndcompareRIN.m b/Time-Series-Analyzer/computeAndcompareRIN.m new file mode 100644 index 0000000..1083007 --- /dev/null +++ b/Time-Series-Analyzer/computeAndcompareRIN.m @@ -0,0 +1,143 @@ +% Script to compute the Relative Intensity Noise of a laser by recording the y-t signal +% by Mathias Neidig in 2012 +% modified for DyLab use by Karthik in 2024 + +% The RIN is defined as +% +% RIN = 10* log10 [Single-sided power spectrum density / (average power)] +% +% and is given in [RIN] = dB/Hz + +clear all +close all + +%% Set the directory where the data is + +dirDCData = 'C:\\Users\\Karthik\\Documents\\GitRepositories\\Calculations\\Time-Series-Analyzer\\Time-Series-Data\\20240915\\After_AOD\\DC_Coupling\\'; +dirACData = 'C:\\Users\\Karthik\\Documents\\GitRepositories\\Calculations\\Time-Series-Analyzer\\Time-Series-Data\\20240915\\After_AOD\\AC_Coupling\\'; + +%% Load the files which contain: - the DC coupled y-t signal to obtain the averaged power +% - the AC coupled y-t signal to obtain the fluctuations +% - the AC coupled y-t signal with the beam blocked to obtain the background fluctuations + +%------------------------------------------------------------------------- % +bgsignal = load( [ dirACData '20240915_Background_AC'] ); % +dcsignal = load( [ dirDCData '20240915_1V_-3Mod_DC_PID_Inactive'] ); % +acsignal_1 = load( [ dirACData '20240915_1V_-3Mod_AC_PID_Inactive'] ); % +acsignal_2 = load( [ dirACData '20240915_1V_-3Mod_AC_PID_Active'] ); % + +label_0 = 'Background'; % +label_1 = 'Power = 1 V, Modulation = -3.0 V, PID Inactive'; % +label_2 = 'Power = 1 V, Modulation = -3.0 V, PID Active'; % +%------------------------------------------------------------------------- % + +%% Read out the important parameters + +dcdata = dcsignal.A; +acdata_1 = acsignal_1.A; +acdata_2 = acsignal_2.A; +bgdata = bgsignal.A; + +N = length(dcdata); % #samples +f_s = 1/dcsignal.Tinterval; % Sample Frequency +delta_f = f_s/N; % step size in frequency domain +delta_t = 1/f_s; % time step + +%% Custom Control Parameters + +% Choose smoothing parameter; has to be odd +%----------------% + span = 21; % +%----------------% + +%% Computes the RIN + +% compute the average power (voltage^2) +average_P = mean(dcdata.*dcdata); + +% compute the power spectrum density FFT(A) x FFT*(A)/N^2 of the source & the bg +psd_src_1 = fft(acdata_1) .* conj(fft(acdata_1))/N^2; +psd_src_2 = fft(acdata_2) .* conj(fft(acdata_2))/N^2; +psd_bg = fft(bgdata) .* conj(fft(bgdata))/N^2; + +% converts the psd to the single-sided psd --> psd is symmetric around zero --> omit +% negative frequencies and put the power into the positive ones --> spsd + +for i = 1 : N/2+1 + if i>1 + spsd_src_1(i) = 2*psd_src_1(i); + spsd_src_2(i) = 2*psd_src_2(i); + spsd_bg(i) = 2*psd_bg(i); + else + spsd_src_1(i) = psd_src_1(i); + spsd_src_2(i) = psd_src_2(i); + spsd_bg(i) = psd_bg(i); + end +end + +% smooths the spsd by doing a moving average +spsd_src_smooth_1 = smooth(spsd_src_1,span,'moving'); +spsd_src_smooth_2 = smooth(spsd_src_2,span,'moving'); +spsd_bg_smooth = smooth(spsd_bg, span,'moving'); + +% calculates the RIN given in dB/Hz; the factor delta_f is needed to convert from dB/bin into dB/Hz +RIN_src_smooth_1 = 10*log10(spsd_src_smooth_1/(average_P*delta_f)); +RIN_src_smooth_2 = 10*log10(spsd_src_smooth_2/(average_P*delta_f)); +RIN_bg_smooth = 10*log10(spsd_bg_smooth /(average_P*delta_f)); + +% creates an array for the frequencies up to half the sampling frequency +f = f_s/2 * linspace(0,1,N/2+1); +f_smooth = smooth(f,span,'moving'); +% +% Calculates the shot noise limit of the used PD given the wavelength of the light source and +% incident average power +PlanckConstant = 6.62607015E-34; +SpeedOfLight = 299792458; +WavelengthOfLaserLight = 1064E-9; +FrequencyOfLaserLight = SpeedOfLight / WavelengthOfLaserLight; +QuantumEfficiencyOfPD = 1; +AverageIncidentPower = 0.001; % (in W) +ShotNoiseLimit = 10*log10((2 * PlanckConstant * FrequencyOfLaserLight) / (QuantumEfficiencyOfPD * AverageIncidentPower)); + +%% Plots the RIN + +% Plots the RIN vs frequency +f_ = clf; +figure(f_); +set(gcf,'Position',[100 100 950 750]) +semilogx(f_smooth, RIN_bg_smooth, LineStyle = "-", Color = [.7 .7 .7]) +hold on +semilogx(f_smooth,RIN_src_smooth_1,'r-') +semilogx(f_smooth,RIN_src_smooth_2,'b-') +ax = gca(f_); +ax.XAxis.FontSize = 14; +ax.YAxis.FontSize = 14; +xlabel('Frequency [Hz]', FontSize=16) +ylabel('RIN [dBc/Hz]', FontSize=16) +xlim([10 5E6]); +ylim([-140 -25]); +title('\bf Intensity noise of Arm 1 of the ODT as measured by an Out-of-loop PD', FontSize=14) +legend(label_0, label_1, label_2, 'Location','NorthWest', FontSize=12); +grid on + +% optional: save the picture without editing wherever you want +%------------------------------------------% +% saveas(f_,'FileName','png'); % +%------------------------------------------% + + + + + + + + + + + + + + + + +