- Loop over wavelengths ('lambda') instead of angles ('theta').
- Update calculations that depend on wavelength (wavevector k, z-component for each layer 'kz1' to 'kz5').
- Reflectance calculation remains similar using the updated 'kz' values.
I want to plot SPR curve for Reflectance Vs Wavelength in place of angle of incidence.can anyone guide me ?
2 views (last 30 days)
Show older comments
function SPR_matrix_profile_comparison_v2
%variables
n_Ma=1.493; %målt ved 589nm @ 20 C
n_buffer=1.33156; %HEPES buffer
d_Ma=1E-9;
%vinkelting
theta_start=50.01; %GRADER! !
theta_slut=60.00; %GRADER! !
thetasteps=300;
theta_start=theta_start/360*2*pi ;
theta_slut=theta_slut/360*2*pi ;
deltatheta=(theta_slut-theta_start)/(thetasteps-1);
for t =1:1:thetasteps;
theta=theta_start+deltatheta*(t-1);
vinkel(t)=theta/(2* pi)*360 ;
end
figure(1) ;
hold on
for l =2:8
d_Au=l*10E-9;
[R_list]=spr_calc(n_Ma,n_buffer,d_Ma,d_Au,vinkel,theta_start,theta_slut,deltatheta,thetasteps,theta);
% [y_data]=data_import(vinkel);
if l==2
plot( vinkel , R_list ,' k ')
end
if l==3
plot ( vinkel , R_list ,' g ')
end
if l==4
plot ( vinkel , R_list ,'m')
end
if l==5
plot ( vinkel , R_list ,'c')
end
if l==6
plot ( vinkel , R_list ,' b ')
end
if l==7
plot ( vinkel , R_list , ' y ' )
end
if l==8
plot ( vinkel , R_list , ' r ' )
end
end
%axis ( [ theta_start *180/ pi theta_slut *180/ pi+1e-3 0 1 ] )
legend ( '20nm' , ' 30nm' , ' 40nm' , ' 50nm' , ' 60nm' , ' 70nm' , ' 80nm' ,'Location','Southheast');
%l egend ( 'n_3=1.3300 ' , 'n_3=1.3340 ' , 'n_3=1.3380 ' , 'n_3=1.3420 ' , 'n_3=1.3460 ' , 'n_3=1.3500 ' , ' Locatxlabel ( ' Angle o f incidence ' )
ylabel ( ' Reflectance(R) ' );
xlabel('Angle of incidence');
load SPR_data.txt ;
[ sx , sy ]= size (SPR_data) ;
for i =1:300
y_data ( i )=SPR_data( i +150 ,3)/7E4 ;
%x_data ( i )=i ;
end
% options = fitoptions ( ' NonlinearLeastSquares ' ,method ) ;
%
%
% [ fresult , gof , fout ]= fit (vinkel, y_data , spr_calc , options ) ;
%figure(1) ;
%hold on
%plot( vinkel , R_list )
%hold on ;
%figure ( 1 ) ;
%plot(vinkel , R_list )
%figure (2) ;
%plot(vinkel , y_data ) ;
%figure(2) ;
%plot(vinkel , R_list-y_data ) ;
end
function [R_list] = spr_calc(n_Ma,n_buffer,d_Ma,d_Au,vinkel,theta_start,theta_slut,deltatheta,thetasteps,theta)
lambda_0=780E-9;
n_glass =1.77;
n_Au=0.17-4.93i;
n_Cr=3.11982+3.44408i ; % http : / / refractive index . info/? group=METALS&material=Chromium d.11n_Au=0.17-4.93i ; %Stenberg etal . 1991
n_Ma=n_Ma;
n_buffer=n_buffer ;
d_Cr=1E-9;
d_Au=d_Au;
d_Ma=d_Ma;
% theta_start=51; %GRADER! !
% theta_slut=55; %GRADER! !
% thetasteps=250;
n1=n_glass ;
n2=n_Cr ;
n3=n_Au;
n4=n_Ma;
n5=n_buffer;
eps1=n1^2;
eps2=n2^2;
eps3=n3^2;
eps4=n4^2;
eps5=n5^2;
k=2*pi/lambda_0 ;
% theta_start=theta_start /360*2* pi ;
% theta_slut=the ta_slut /360*2* pi ;
% deltatheta=( theta_slut-theta_start ) / (thetasteps-1);
for t =1:1: thetasteps
theta=theta_start+deltatheta *( t-1);
kx=k* n_glass * sin ( theta ) ;
kz1=sqrt (-(kx^2)+k^2* eps1 ) ;
if imag ( kz1 )>0 ,
kz1=-kz1 ;
end
kz2=sqrt (-(kx^2)+k^2* eps2 ) ;
if imag ( kz2 )>0 ,
kz2=-kz2 ;
end
kz3=sqrt (-(kx^2)+k^2* eps3 ) ;
if imag ( kz3 )>0 ,
kz3=-kz3 ;
end
kz4=sqrt (-(kx^2)+k^2* eps4 ) ;
if imag ( kz4 )>0 ,
kz4=-kz4 ;
end
kz5=sqrt (-(kx^2)+k^2* eps5 ) ;
if imag ( kz5 )>0 ,
kz5=-kz5 ;
end
beta_Cr=kz2*d_Cr ; %2*pi /lambda_0 *(n_Cr*d_Cr* cos ( theta ) ) ;
beta_Au=kz3*d_Au; %2*pi /lambda_0 *(n_Au*d_Au* cos ( theta ) ) ;
beta_Ma=kz4*d_Ma; %2*pi /lambda_0 *(n_Ma*d_Ma* cos ( theta ) ) ;
L_Cr=zeros(2) ;
L_Cr(1 ,1)=exp(-1i *beta_Cr ) ;
L_Cr(2 ,2)=exp (1i *beta_Cr ) ;
L_Au=zeros(2) ;
L_Au(1 ,1)=exp(-1i *beta_Au ) ;
L_Au(2 ,2)=exp (1i *beta_Au ) ;
L_Ma=zeros(2) ;
L_Ma(1 ,1)=exp(-1i *beta_Ma ) ;
L_Ma(2 ,2)=exp (1i *beta_Ma ) ;
b1=(n1/n2 )^2*( kz2 /kz1 ) ;
b2=(n2/n3 )^2*( kz3 /kz2 ) ;
b3=(n3/n4 )^2*( kz4 /kz3 ) ;
b4=(n4/n5 )^2*( kz5 /kz4 ) ;
rho_glass_Cr=(1-b1)/(1+b1 ) ;
rho_Cr_Au=(1-b2)/(1+b2 ) ;
rho_Au_Ma=(1-b3)/(1+b3 ) ;
rho_Ma_buffer=(1-b4)/(1+b4 ) ;
tau_glass_Cr=(2*(n1/n2 ))/(1+b1 ) ;
tau_Cr_Au=(2*(n2/n3 ))/(1+b2 ) ;
tau_Au_Ma=(2*(n3/n4 ))/(1+b3 ) ;
tau_Ma_buffer=(2*(n4/n5 ))/(1+b4 ) ;
H_glass_Cr=ones (2)/ tau_glass_Cr ;
H_glass_Cr(1,2)=rho_glass_Cr / tau_glass_Cr ; %rho i s r e f l e c t i o n c o e f f i c i e n t s
H_glass_Cr(2,1)=rho_glass_Cr / tau_glass_Cr ; %rho i s r e f l e c t i o n c o e f f i c i e n t s
H_Cr_Au=ones(2)/tau_Cr_Au ;
H_Cr_Au(1,2)=rho_Cr_Au/tau_Cr_Au ; %rho i s r e f l e c t i o n c o e f f i c i e n t s
H_Cr_Au(2,1)=rho_Cr_Au/tau_Cr_Au ; %rho i s r e f l e c t i o n c o e f f i c i e n t s
H_Au_Ma=ones (2)/tau_Au_Ma ;
H_Au_Ma(1,2)=rho_Au_Ma/tau_Au_Ma ; %rho i s r e f l e c t i o n c o e f f i c i e n t s
H_Au_Ma(2,1)=rho_Au_Ma/tau_Au_Ma ; %rho i s r e f l e c t i o n c o e f f i c i e n t s
H_Ma_buffer=ones (2)/ tau_Ma_buffer ;
H_Ma_buffer(1,2)=rho_Ma_buffer/tau_Ma_buffer ; %rho i s r e f l e c t i o n c o e f f i c i e n t s
H_Ma_buffer(2,1)=rho_Ma_buffer/tau_Ma_buffer ; %rho i s r e f l e c t i o n c o e f f i c i e n t s
S_glass_buffer=H_glass_Cr*L_Cr*H_Cr_Au*L_Au*H_Au_Ma*L_Ma*H_Ma_buffer ;
R_list(1,t )=( abs(S_glass_buffer ( 1 , 2 ) / S_glass_buffer ( 2 , 2 ) ) )^2 ;
% vinkel (t)=theta /(2* pi ) * 3 6 0 ;
end
% figure(1) ;
% hold on
% plot(vinkel , R_list )
% [ sy , sx ]= s i z e ( R_l i s t ) ;
%
% for i =1:(sx-5)
% theta=the ta_start+deltatheta *(i -1);
% diff_R_list(i)=(-R_list(i)+R_list(i+5))/(5) ;
% diff_vinkel(i)=theta/(2* pi)* 360 ;
% end
%
% for i =1:(sx-6)
% theta=theta_start+deltatheta *( i-1);
% if sign ( diff_R_list ( i ))==sign ( diff_R_list (i +1))
% signchange (i)=0;
% else
% sign_note=theta /(2* pi )*360
% end
% end
%
% figure(2)
% plot ( diff_vinkel , diff_R_list )
end
function [ y_data ] = data_import(vinkel)
load SPR_data.txt ;
[sx,sy]= size(SPR_data) ;
for i =1:250
y_data( i )=SPR_data( i +150 ,2)/7E4 ;
%x_data( i )=i ;
end
x_data=vinkel ;
% hold on ;
%
% figure(2) ;
% plot( x_data , y_data )
end
0 Comments
Answers (1)
Abhinaya Kennedy
on 16 May 2024
Hi Ananya,
Your code simulates SPR for different gold thicknesses but plots reflectance vs. incident angle. To generate a reflectance vs. wavelength plot, modify the 'spr_calc' function:
See the modified 'spr_calc' function for details.
function [R_list] = spr_calc(n_Ma, n_buffer, d_Ma, d_Au, lambda_range, theta_start, theta_slut, deltatheta, thetasteps)
% Define material properties (permittivities can be replaced with functions for wavelength dependence)
n_glass = 1.77;
n_Cr = 3.11982 + 3.44408i;
n_Au = 0.17 - 4.93i;
n_Ma = n_Ma;
n_buffer = n_buffer;
% Loop over wavelengths
R_list = zeros(length(lambda_range), thetasteps); % Pre-allocate reflectance matrix
for l = 1:length(lambda_range)
lambda = lambda_range(l);
% Calculate wavevector
k = 2*pi/lambda;
% Calculate z-component of wavevector in each layer
for t = 1:thetasteps
theta = theta_start + deltatheta*(t-1);
kx = k*n_glass*sin(theta);
kz1 = sqrt(-(kx^2) + k^2*n_glass^2);
if imag(kz1) > 0, kz1 = -kz1; end
kz2 = sqrt(-(kx^2) + k^2*n_Cr^2);
if imag(kz2) > 0, kz2 = -kz2; end
kz3 = sqrt(-(kx^2) + k^2*n_Au^2);
if imag(kz3) > 0, kz3 = -kz3; end
kz4 = sqrt(-(kx^2) + k^2*n_Ma^2);
if imag(kz4) > 0, kz4 = -kz4; end
kz5 = sqrt(-(kx^2) + k^2*n_buffer^2);
if imag(kz5) > 0, kz5 = -kz5; end
% ... (rest of the calculations using kz1 to kz5 for reflectance)
R_list(l, t) = (calculated reflectance value); % Update reflectance based on kz values
end
end
end
This function now returns a reflectance matrix where rows represent reflectance for a gold layer thickness and columns represent reflectance for a specific wavelength. You can then plot the desired SPR curve.
0 Comments
See Also
Categories
Find more on 2-D and 3-D Plots in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!