%《MATLAB及其在理工課程中的應用指南》第六章例6－7－1程序
% 單色光的雙縫干涉
% 西安電子科技大學出版社出版  陳懷琛編著 1999年10月初版，2004年11月第二版
%
clear; 
Lambda = input('輸入光的波長(單位為 nm): (書上取500) ');
Lambda = Lambda * 1e-9;                 % 將nm換 變為 m
d = input('輸入兩個縫的間距 (單位為 mm): (書上取2) ');
d = d * 0.001;  				        % 將mm 變換為 m
Z = input('輸入縫到屏的距離 (單位為 m): (書上取1) ');
yMax = 5*Lambda*Z/d; xs = yMax;         % 設定圖案的y,x 向范圍
Ny=101;ys = linspace(-yMax,yMax,Ny);    % y方向分成101點
for i=1:Ny			                    % 對屏上全部點進行循環計算
    % 計算第一和第二個光源到屏上各點的距離
    L1 = sqrt((ys(i)-d/2).^2 + Z^2 );
    L2 = sqrt((ys(i)+d/2).^2 + Z^2 );
    Phi = 2*pi*(L2-L1)/Lambda;          % 從距離差計算相位差
    B(i,:) = 4*cos(Phi/2).^2;           % 計算該點光強（設兩束光強相同)
    % 若考慮光譜的非單色性,把前兩句改為后四句
    %Nl=11; dL=linspace(-0.1,0.1,Nl);	% 設光譜相對寬度正負10%,
    %Lambda1=Lambda*(1+dL');			% 分11根譜線,波長為一個數組
    %Phi1 = 2*pi*(L2-L1)./Lambda1;      % 從距離差計算各波長的相位差
    %B(i,:) = sum(4*cos(Phi1/2).^2)/Nl; % 疊加各波長影響,計算光強
end
% 在屏上畫出圖象
% clf; figure(gcf);                     % 清圖形窗，將它移到前面
NCLevels = 255;                         % 確定用的灰度等級
% 定標：使最大光強(4.0)對應于最大灰度級（白色）
Br = (B/4.0) * NCLevels;
subplot(1,4,1),image(xs,ys,Br);         % 畫圖象
colormap(gray(NCLevels));               % 用灰度級顏色圖
subplot(1,4,2), plot(B(:),ys)           % 畫出沿y 向的光強變化曲線