Updated to latest versions of this codeHEADmaster

4 files changed, 65 insertions, 47 deletions

diff --git a/getPeak.m b/getPeak.m
index 088a9ca..b074c73 100644
--- a/getPeak.m
+++ b/getPeak.m
@@ -1,4 +1,4 @@
-function [fwhm, absorption_min] = getPeak(detunings, absorptions)
+function [fwhm, absorption_min, contrast_left, contrast_right, features] = getPeak(detunings, absorptions)
 
     eit_data = [detunings, absorptions];
     
@@ -27,16 +27,24 @@ function [fwhm, absorption_min] = getPeak(detunings, absorptions)
     absorption_mid = (absorption_max + absorption_min) / 2;
     half_width_right = interp1(eit_data(right_range,2), eit_data(right_range,1), absorption_mid) - eit_data(peak_center(1),1);
     half_width_left = eit_data(peak_center(1),1) - interp1(eit_data(left_range,2), eit_data(left_range,1), absorption_mid);
-    fwhm = half_width_right + half_width_left;
-    
     
+    fwhm = half_width_right + half_width_left;
+    linewidth = [eit_data(peak_center(1),1)-half_width_left, absorption_mid; eit_data(peak_center(1),1)+half_width_right, absorption_mid];
+    left_max = [detunings(left_max(1)), left_max(2)];
+    right_max = [detunings(right_max(1)), right_max(2)];
+    peak_center = [detunings(peak_center(1)), peak_center(2)];
+    features = [left_max; peak_center; right_max; linewidth];
+    contrast_left = (left_max(2) - peak_center(2))/left_max(2);
+    contrast_right = (right_max(2) - peak_center(2))/right_max(2);
+
     %Plot extrema over existing plot
-    % plot(eit_data(:,1), eit_data(:,2));
-    % hold all
-    % grid on
-    % plot(detuning_1(left_max(1)), left_max(2), 'r*')
-    % plot(detuning_1(right_max(1)), right_max(2), 'r*')
-    % plot(detuning_1(peak_center(1)), peak_center(2), 'g*')
-    
+    % features = figure('visible', 'off');
+    % plot(
+    %     eit_data(:,1), eit_data(:,2),
+    %     detunings(left_max(1)), left_max(2), 'r*',
+    %     detunings(right_max(1)), right_max(2), 'r*',
+    %     detunings(peak_center(1)), peak_center(2), 'g*',
+    %     linewidth(:,1), linewidth(:,2)
+    % );
 
 end
diff --git a/loadSimulation.m b/loadSimulation.m
index e5f0be3..2626cdb 100644
--- a/loadSimulation.m
+++ b/loadSimulation.m
@@ -2,6 +2,10 @@ function [width, minAbsorption, detunings, absorptions] = loadSimulation(dataFil
 	data = load(dataFile);
 	detunings = data(:,1);
 	absorptions = data(:,2);
-	[width, minAbsorption] = getPeak(detunings, absorptions);
-	plot(detunings, absorptions);
+	[width, minAbsorption, transmission, features] = getPeak(detunings, absorptions);
+	figure
+	hold all
+	plot(detunings, absorptions)
+	plot(features(:,1), features(:,2), 'r*')
+
 end
\ No newline at end of file
diff --git a/loadSimulations.m b/loadSimulations.m
index 046e336..306fe7e 100644
--- a/loadSimulations.m
+++ b/loadSimulations.m
@@ -1,27 +1,30 @@
-function [data] = loadSimulations(decay_bc, dephase_bc)
-	dataFiles = dir(['data/', decay_bc, '*_*_*_', dephase_bc, '*.dat']);
+function [data] = loadSimulations(directory, decay_bc, detuning, dephase_bc)
+	dataFiles = dir([directory, decay_bc, '*_*_', detuning, '*_', dephase_bc, '*.dat']);
 	widths = [];
 	minAbsorptions = [];
+	contrasts = [];
 	drives = [];
 	deltas = [];
 	for i = 1:length(dataFiles)
-		dataFile = ['data/', dataFiles(i).name];
+		dataFile = [directory, dataFiles(i).name];
 		data = load(dataFile);
 		detunings = data(:,1);
 		absorptions = data(:,2);
-		[width, absorption] = getPeak(detunings, absorptions);
+		[width, absorption, contrast_left, contrast_right] = getPeak(detunings, absorptions);
+		contrast = [contrast_left, contrast_right];
 		widths = [widths; width];
 		minAbsorptions = [minAbsorptions; absorption];
+		contrasts = [contrasts; contrast];
 		sections = strfind(dataFile, '_');
 		drive = dataFile(sections(1)+1:sections(2)-1);
 		drives = [drives; str2num(drive)];
 		delta = dataFile(sections(2)+1:sections(3)-1);
 		deltas = [deltas; str2num(delta)];
 	end
-	data = [drives, widths, minAbsorptions];
+	data = [drives, widths, minAbsorptions, contrasts];
 	data = sortrows(data, 1);
-	drives = data(:,1); widths = data(:,2); minAbsorptions = data(:,3);
-	loglog(drives.^2, widths);
+	drives = data(:,1); widths = data(:,2); minAbsorptions = data(:,3); contrasts = data(:,4:5);
+	loglog(drives, widths);
 	hold all
-	loglog(drives.^2, minAbsorptions*10^5);
+	loglog(drives, minAbsorptions.*10^5);
 end
\ No newline at end of file
diff --git a/simulate.m b/simulate.m
index ebf882f..ba34bd5 100644
--- a/simulate.m
+++ b/simulate.m
@@ -7,42 +7,45 @@ setenv('LD_LIBRARY_PATH', libPath);
 system('xmds2 eit.xmds');
 setenv('LD_LIBRARY_PATH', matlabPath);
 
-decay_bc = 0.001;
+decay_bc = 0.00001;
 decay_ab = 6;
 
 widths = [];
 absorptions = [];
-drives = logspace(-2, -3, 10);
-drives = linspace(0.003298,0.003319,20);
-detunings1 = [0:2:0];
-dephase_bc = 0.05;
+drives = logspace(0, -2, 100);
+drives = drives(1:1);
+%drives = linspace(0.0033154997,0.00331549970606061,100);
+detunings1 = [-15];
+dephase_bc = [0];
 
 for drive = drives
     width_row = [];
     absorption_row = [];
-    expected_width = (drive^2/decay_ab + dephase_bc + decay_bc);
-    domain_max = expected_width*10;
-    domain_min = -1*domain_max;
-    final_time = 1/expected_width*100;
     for detuning1 = detunings1
-    
-        run_eit = sprintf('time ./eit --drive_arg=%f --decay_arg=%f --decay_bc_arg=%f --delta1_arg=%f --domain_min=%f --domain_max=%f --dephase_bc_arg=%f --final_time=%f', drive, decay_ab, decay_bc, detuning1, domain_min, domain_max, dephase_bc, final_time);
-        
-        setenv('LD_LIBRARY_PATH', libPath);
-        system(run_eit);
-        system('xsil2graphics2 -m eit.xsil');
-        setenv('LD_LIBRARY_PATH', matlabPath);
-        %[drive, detuning1] % Use for debugging
-        run('eit.m');
-        detunings2 = detuning_1';
-        PabI = PabI_1(end,:)';
-        [fwhm, absorption_min] = getPeak(detunings2, PabI);
-        
-        width_row = [width_row, fwhm];
-        absorption_row = [absorption_row, absorption_min];
-        eit_plot = [detunings2, PabI];
-        csv_name = sprintf('data/%f_%f_%f_%f.dat', decay_bc, drive, detuning1, dephase_bc);
-        csvwrite(csv_name, eit_plot);
+        for dephasing = dephase_bc
+
+            expected_width = (drive^2/decay_ab + dephasing + decay_bc);
+            domain_max = expected_width*10;
+            domain_min = -1*domain_max;
+            final_time = 1/expected_width*100;
+            run_eit = sprintf('time ./eit --drive_arg=%0.20f --decay_arg=%f --decay_bc_arg=%f --delta1_arg=%f --domain_min=%f --domain_max=%f --dephase_bc_arg=%f --final_time=%f', drive, decay_ab, decay_bc, detuning1, domain_min, domain_max, dephasing, final_time);
+            
+            setenv('LD_LIBRARY_PATH', libPath);
+            system(run_eit);
+            system('xsil2graphics2 -m eit.xsil');
+            setenv('LD_LIBRARY_PATH', matlabPath);
+            %[drive, detuning1] % Use for debugging
+            run('eit.m');
+            detunings2 = detuning_1';
+            PabI = PabI_1(end,:)';
+            [fwhm, absorption_min] = getPeak(detunings2, PabI);
+            
+            width_row = [width_row, fwhm];
+            absorption_row = [absorption_row, absorption_min];
+            eit_plot = [detunings2, PabI];
+            csv_name = sprintf('data/full/%0.10f_%0.20f_%0.10f_%0.10f.dat', decay_bc, drive, detuning1, dephasing);
+            csvwrite(csv_name, eit_plot);
+        end
     end
     widths = [widths; width_row];
     absorptions = [absorptions; absorption_row];