MATLAB：杨氏双缝干涉仿真程序

念着倒啥为

已于 2023-12-21 22:11:06 修改

阅读量787

点赞数 11

文章标签： matlab

于 2023-12-18 19:56:52 首次发布

本文链接：https://blog.csdn.net/frlove/article/details/135069772

版权

原理

其中光程差计算：

$\Delta L=\frac{x}{D}\cdot d$

光强分布公式：

$I=4I_{0}cos^{^{2}}(\frac{\pi\cdot \Delta L}{\lambda } )$

程序代码

clc
close(findall(0, 'type', 'figure'))
clear;
interference()

function interference()
%% 初始参数(可更改)
stepd=0.01;
stepD=3;
mind=0.1;  maxd=0.5;
minD=500;    maxD=1000;
X=10;      % 显示范围[-10 10]
default_d=mind+round((maxd-mind)/2,1);  % 双缝间距d的默认值
default_D=minD+round((maxD-minD)/2,1);   % 屏间距D的默认值
default_lam=650; % 默认波长
%% 界面
Fig=uifigure('Position',[100,500,800,400],'MenuBar','none');
Fig.Color='white';
Pan1=uipanel('Parent',Fig,'Units','normalized','Position',[0.05,0.4,0.9,0.6]);
Pan2=uipanel('Parent',Fig,'Units','normalized','Position',[0.05,0.03,0.9,0.35]);
Axes1=uiaxes(Pan1,'Position',[8,15,700,180]);
uilabel(Pan1,'Position',[320 200 150 40],'Text','双缝干涉','FontSize',24,'FontWeight','bold');
%% 数值编辑框（步长）
edtd=uieditfield(Pan2,'numeric','Position',[230 80 50 30],'Value',stepd,'HorizontalAlignment','center', ...
    'ValueChangedFcn',@(edtd,event) dValChanged(event));
edtD=uieditfield(Pan2,'numeric','Position',[230 15 50 30],'Value',stepD,'HorizontalAlignment','center', ...
    'ValueChangedFcn',@(edtD,event) DValChanged(event));

%% 滑块d
sld_d = uislider(Pan2,'Position',[300 100 300 3],'Value',default_d,'Limits',[mind maxd]);
Spilad= uispinner(Pan2,'Enable','on','Editable','on','Position',[620 80 70 30],'Value',default_d, ...
    'Step',edtd.Value,'Limits',[mind maxd]);
uilabel(Pan2,'Position',[300 110 150 20],'Text','缝间距d(单位mm)','FontSize',14);

%% 滑块D
sld_D = uislider(Pan2,'Position',[300 30 300 3],'Value',default_D,'Limits',[minD maxD]);
SpilaD= uispinner(Pan2,'Enable','on','Editable','on','Position',[620 10 70 30],'Value',default_D, ...
    'Step',edtD.Value,'Limits',[minD maxD]);
uilabel(Pan2,'Position',[300 40 150 20],'Text','屏间距D(单位mm)','FontSize',14);

%% 按钮
uibutton(Pan2, 'Text', '设置','Position',[60,20, 100, 30],'FontWeight','bold','FontSize',16, ...
    'ButtonPushedFcn', @(btn,event) Set());

%% 波长微调器
Spilam= uispinner(Pan2,'Enable','on','Editable','on','Position',[60 65 100 30], ...
    'Value',default_lam,'FontSize',16,'FontWeight','bold','HorizontalAlignment','center','ValueDisplayFormat','%.0f nm');
uilabel(Pan2,'Position',[65 100 150 20],'Text','波长','FontSize',16);


%% 回调函数定义
Spilam.ValueChangingFcn=@(Spilam,event) lamspinnerChanging(event);
sld_d.ValueChangedFcn=@(sld_d,event) ddateGauge(sld_d);
sld_D.ValueChangedFcn=@(sld_D,event) DdateGauge(sld_D);
Spilad.ValueChangingFcn=@(Spilad,event) dspinnerChanging(event);
SpilaD.ValueChangingFcn=@(SpilaD,event) DspinnerChanging(event);
%% 按钮
Setc=0;
 function Set()         
        if Setc==0
        opt=[mind,maxd,minD,maxD,X];
        Setc=inputc(Fig,opt);
        else
            Setc=inputc(Fig,Setc);
        end

        if Setc==0
            
        else   
            set(sld_d,'Limits',[Setc(1) Setc(2)])      
            set(Spilad,'Limits',[Setc(1) Setc(2)]) 
            set(sld_d,'Value',Setc(1)+round((Setc(2)-Setc(1))/2,1)) 
            set(Spilad,'Value',Setc(1)+round((Setc(2)-Setc(1))/2,1)) 

            set(sld_D,'Limits',[Setc(3) Setc(4)])      
            set(SpilaD,'Limits',[Setc(3) Setc(4)])
            set(sld_D,'Value',Setc(3)+round((Setc(4)-Setc(3))/2,1)) 
            set(SpilaD,'Value',Setc(3)+round((Setc(4)-Setc(3))/2,1)) 
            X=Setc(5);
            showimage();
        end
    end
%% 双峰干涉
 lambda=Spilam.Value;     % 波长（单位：nm）  
showimage();
%% 回调函数   
    function showimage()
        I0=2;  % 光强                
        d=sld_d.Value;     % 双缝间距（单位：mm）
        D=sld_D.Value;      % 双缝到显示间距（单位：mm）
        x=-X:0.01:X;   % 显示范围（单位：mm）
        lam=lambda*10^(-6); 
        %I=@(x)2*I0*(1+cos((2*pi.*x*d)/(lam*D))); 
        I=@(x)4*I0*cos((pi.*x*d)/(lam*D)).^2;
        kx=(D/d)*lam;       
        k=floor(X/kx);
        xk=(-k:1:k)*kx;
        Ix=I(x);
        Norm_I=rescale(Ix);
        [r,g,b] = lam2rgb(lambda);  
        num=length(Norm_I);
        im=ones(round(num/10),num,3);
        im(:,:,1)=r;
        im(:,:,2)=g;
        im(:,:,3)=b;
        for j=1:num
            im(:,j,:)= im(:,j,:)*Norm_I(j);
        end
        image(im,'Parent',Axes1,'XData',[-X,X],'YData',[0,1]);
        axis(Axes1,[-X-0.1 X+0.1 0 1])
        set(Axes1,'xtick',round(xk,3))
        set(Axes1,'ytick',[],'ycolor','w')
    end
    function ddateGauge(event)        
        set(Spilad,'Value',round(event.Value,3))
        showimage();
    end
    function DdateGauge(event)        
        set(SpilaD,'Value',round(event.Value,3))
        showimage();
    end

    function lamspinnerChanging(event)
        if event.Value<380
            uialert(Fig,'请输入380-780nm范围的波长','警告');
            pause(0.1)
            set(Spilam,'Value',380)

        elseif event.Value>780
            uialert(Fig,'请输入380-780nm范围的波长','警告');
            pause(0.1)
            set(Spilam,'Value',780)
        else
            lambda=event.Value;
            showimage();
        end
    end
    function dspinnerChanging(event)
        set(sld_d,'Value',event.Value)
        
        showimage();
    end
    function DspinnerChanging(event)
        set(sld_D,'Value',event.Value)
       
        showimage();
    end
    function dValChanged(event)
        set(Spilad,'Step',event.Value)
    end
    function DValChanged(event)
        set(SpilaD,'Step',event.Value)
    end


%% 波长转化RGB
    function [r,g,b] = lam2rgb(lambda)

        if  (lambda >= 380.0) && (lambda < 440.0)
            r = -1.0 * (lambda - 440.0) / (440.0 - 380.0);
            g = 0.0;
            b = 1.0;
        elseif (lambda >= 440.0) && (lambda < 490.0)
            r = 0.0;
            g = (lambda - 440.0) / (490.0 - 440.0);
            b = 1.0;
        elseif (lambda >= 490.0) && (lambda < 510.0)
            r = 0.0;
            g = 1.0;
            b = -1.0 * (lambda - 510.0) / (510.0 - 490.0);
        elseif (lambda >= 510.0) && (lambda < 580.0)
            r = (lambda - 510.0) / (580.0 - 510.0);
            g = 1.0;
            b = 0.0;
        elseif (lambda >= 580.0) && (lambda < 645.0)
            r = 1.0;
            g = -1.0 * (lambda - 645.0) / (645.0 - 580.0);
            b = 0.0;
        elseif (lambda >= 645.0) && (lambda <= 780.0)
            r = 1.0;
            g = 0.0;
            b = 0.0;
        else
            r = 0.0;
            g = 0.0;
            b = 0.0;
        end

        % 在可见光谱的边缘处强度较低。
        if  (lambda >= 380.0) && (lambda < 420.0)
            attenuation = 0.3 + 0.7 * (lambda - 380) / (420 - 380);
            r = r * attenuation;
            g = 0.0;
            b = 1.0 * attenuation;
        elseif(lambda >= 701.0) && (lambda <= 780.0)
            attenuation = 0.30 + 0.70 * (780.0 - lambda) / (780.0 - 700.0);
            r = r * attenuation;
            g = 0.0;
            b = 0.0;
        end
    end
end

function [c]=inputc(fig,daf)
prompt = {'输入双缝间距d下限:','输入双缝间距d上限:','屏间距D下限:','屏间距D上限:','显示范围X'};
dlgtitle = '输入参数';
dims = [1 35];
definput{5}=[];
if isempty(daf)
    definput={'0.1','0.5','500','1000','10'};
else
    for n=1:5
        definput{n}=num2str(daf(n));
    end
end

answer = inputdlg(prompt,dlgtitle,dims,definput);

if isempty(answer)
   c=str2double(definput);
else
    dmin=str2double(answer{1});
    dmax=str2double(answer{2});
    Dmin=str2double(answer{3});
    Dmax=str2double(answer{4});
    X=str2double(answer{5});
    if dmin>=dmax || Dmin>=Dmax
        uialert(fig,'输入最小值大于等于最大值','警告');
        if isempty(daf)
            c=0;
        else
            c=daf;
        end
    else
        c=[dmin,dmax,Dmin,Dmax,X];
    end

end
end