How to find the peaks?

How to find the peak values of the curve. I try, but not getting result.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This is a program for 1-D Photonic crystal
% For Photonic band structure calculation and Transmission Spectrum
% This program uses Transfer Matrix Method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clc
clear all
close all
tic
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Initialization of parameters
dB=7778e-9;                         % Thickness of Second Layer
dA=5303e-9;                     % Thickness of First Layer

nA=3.3;                             % Refractive index of First Layer
nB=2.25;                            % Refraive index of Second Layer
n0=1;                               % Refractive index of Incident and transmission medium
c=3e8;                              % Velocity of Light (m/s)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Frequency Range
frequncy=0.1:0.001:50;                    % Frequncy in THz
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for loop=1:length(frequncy)
    % angular frequency
    w=2*pi*frequncy*1e12;
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    d=dA+dB;    DA=((w(loop))/c)*dA*nA;    DB=((w(loop))/c)*dB*nB;
    m11=cos(DA);   m12=-1i*sin(DA)/nA;   m21=-1i*nA*sin(DA);   m22=cos(DA);   mA=[m11 m12; m21 m22];
    l11=cos(DB);   l12=-1i*sin(DB)/nB;   l21=-1i*nB*sin(DB);   l22=cos(DB);   mB=[l11 l12; l21 l22];
    Da=((w(loop))/c)*(dA/2)*nA;    Db=((w(loop))/c)*(dB/2)*nB;
    ma11=cos(Da);   ma12=-1i*sin(Da)/nA;   ma21=-1i*nA*sin(Da);   ma22=cos(Da);   ma=[ma11 ma12; ma21 ma22];
    lb11=cos(Db);   lb12=-1i*sin(Db)/nB;   lb21=-1i*nB*sin(Db);   lb22=cos(Db);   mb=[lb11 lb12; lb21 lb22];

    % Unit cell
    m10=(mA*mB);
    
    % Whole 1-D PC structure

    M10=(mb*mA*mb)^2*(ma*mB*ma)^2;
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % Reflection and Transmission

    T10(loop)=(2*n0/(n0*M10(1,1)+n0*n0*M10(1,2)+M10(2,1)+n0*M10(2,2)));Transmission10(loop)=(T10(loop)*conj(T10(loop)));
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    g=double(Transmission10(:));
    [t_max,loc]=findpeaks(g);
  
end

Error using findpeaks>parse_inputs
Data set must contain at least 3 samples.

Error in findpeaks (line 135)
= parse_inputs(isInMATLAB,Yin,varargin{:});

plot(frequncy,Transmission10(:),'LineWidth',1);hold on;
plot(frequncy,double(Transmission10),'or');hold off
ylim([0,1])
toc

[akshatsood]

Hi GULZAR,

I understand that you are facing an error while determining the peak values of a curve. I investigated your code and got an error message which points towards invalid input being supplied to the findpeaks function. As per the documentation page, the input data must be real and needs to have at least 3 elements to determine the peak values of a curve. In your case, you have supplied g to the findpeaks function during the first iteration of the loop which is a single double value causing the error you observe on running the code.

To address this, you can shift the following two statements after the for loop to compute the peak values of the curve.

g=double(Transmission10(:));
[t_max,loc]=findpeaks(g);

I hope this helps.

Answer 1

[Walter Roberson]

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This is a program for 1-D Photonic crystal
% For Photonic band structure calculation and Transmission Spectrum
% This program uses Transfer Matrix Method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clc
clear all
close all
tic
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Initialization of parameters
dB=7778e-9;                         % Thickness of Second Layer
dA=5303e-9;                     % Thickness of First Layer
nA=3.3;                             % Refractive index of First Layer
nB=2.25;                            % Refraive index of Second Layer
n0=1;                               % Refractive index of Incident and transmission medium
c=3e8;                              % Velocity of Light (m/s)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Frequency Range
frequncy=0.1:0.001:50;                    % Frequncy in THz

There is no assignment to Transmission10 up to this point. The variable does not exist before the loop.

for loop=1:length(frequncy)
    % angular frequency
    w=2*pi*frequncy*1e12;
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    d=dA+dB;    DA=((w(loop))/c)*dA*nA;    DB=((w(loop))/c)*dB*nB;
    m11=cos(DA);   m12=-1i*sin(DA)/nA;   m21=-1i*nA*sin(DA);   m22=cos(DA);   mA=[m11 m12; m21 m22];
    l11=cos(DB);   l12=-1i*sin(DB)/nB;   l21=-1i*nB*sin(DB);   l22=cos(DB);   mB=[l11 l12; l21 l22];
    Da=((w(loop))/c)*(dA/2)*nA;    Db=((w(loop))/c)*(dB/2)*nB;
    ma11=cos(Da);   ma12=-1i*sin(Da)/nA;   ma21=-1i*nA*sin(Da);   ma22=cos(Da);   ma=[ma11 ma12; ma21 ma22];
    lb11=cos(Db);   lb12=-1i*sin(Db)/nB;   lb21=-1i*nB*sin(Db);   lb22=cos(Db);   mb=[lb11 lb12; lb21 lb22];
    % Unit cell
    m10=(mA*mB);
    
    % Whole 1-D PC structure
    M10=(mb*mA*mb)^2*(ma*mB*ma)^2;
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % Reflection and Transmission

You start the for loop, loop control variable named loop starts out as 1. You do not assign to Transmission10 between the beginning of the loop and this point.

        T10(loop)=(2*n0/(n0*M10(1,1)+n0*n0*M10(1,2)+M10(2,1)+n0*M10(2,2)));Transmission10(loop)=(T10(loop)*conj(T10(loop)));

You assign to Transmission10(loop) there. The first time through the loop, loop is 1, so you are assigning to Transmission10(1) . The variable did not exist before this point, so it is created here, and indexed by 1 so the size of Transmission10 after this statement in the first iteration of the loop is length 1.

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    g=double(Transmission10(:));

This is still within the loop, so Transmission10 is still length 1, so all one element of it is assigned to g so when loop is 1, g is going to be a scalar.

        [t_max,loc]=findpeaks(g);

and you try to findpeaks() on the scalar g . But findpeaks has a limitation that it refuses to find peaks if the signal is empty or length 1 or length 2.

Are you certain that you want to be doing the findpeaks inside the for loop loop ?