Histogram for a loop

Hi,

I am trying to simulate a chromatogram for 100 molecules using following code and the output looks very strange. What I am expecting is a Gaussian type distribution.

TP=10000; %number of total theoritical plates
N=100; %number of molecules
Pa=0.11; %partition coefficient of A (mobile/stationary)
Pb=0.12; %partition coefficient of B (mobile/stationary)

%This part simulates the retention time for A
tA=zeros(N,1);
for i=1:N %for loop for N number of molecules
NTPa=0; %# of theoritical plates before equilibration
while NTPa<TP; %number of TPs for A should never exceed total number of TPs.
if rand(1,1)>Pa; %partition coefficient of A
NTPa=1+NTPa;
tA(i)=(1+NTPa)*0.01; %tA=retention time of A 
end
end
end

%This part simulates the retention time for B
tB=zeros(N,1);
for i=1:N %for loop for N number of molecules
NTPb=0; %# of theoritical plates before equilibration
while NTPb<TP;
if rand(1,1)>Pb; %partition coefficient of B
NTPb=1+NTPb;
tB(i)=(1+NTPb)*0.01; %tB=retention time of B
end
end
end
subplot(2,1,1);hist(tA/60); %conversion to minutes
xlabel('retention time (min)')
ylabel('Signal (arbitary units)')
title('chromatogram A')
subplot(2,1,2);hist(tB/60); %conversion to minutes
xlabel('retention time (min)')
ylabel('Signal (arbitary units)')
title('chromatogram B')

[Geoff Hayes]

Chanaka - why do you expect to see a Gaussian type distribution? What do you see instead?

[Chanaka Navarathna]

NTPa and NTPb are not going to be a single number right?

I am getting only a single bar on the histogram.

When I do it this way it gives me a Gaussian distribution, but answer is incorrect in that case. I am expecting answers above 1.6 mins based on some theortical calculations I have done to verify.

TP=10000; %number of total theoritical plates
N=1000; %number of molecules
Pa=0.1; %partition coefficient of A (mobile/stationary)
Pb=0.2; %partition coefficient of B (mobile/stationary)

%This part simulates the retention time for A
tA=zeros(N,1);
tB=zeros(N,1);
for i=1:N %for loop for N number of molecules
NTPa=0; %number of theoritical plates of A before equilibration
NTPb=0; %number of theoritical plates of B before equilibration
while NTPa+NTPb<TP; %number of TPs for A+B should never exceed total number of TPs.
if rand(1,1)>Pa; %condition to travel A in the mobile phase
    NTPa=1+NTPa;
if rand(1,1)>Pb; %condition to travel B in the mobile phase
    NTPb=1+NTPb;
tA(i)=(1+NTPa)*0.01+10; %tA=retention time of A
tB(i)=(1+NTPb)*0.01+10; %tB=retention time of B 
end
end
end
end

subplot(2,1,1);hist(tA/60); %conversion to minutes
xlabel('retention time (min)')
ylabel('Signal (arbitary units)')
title('chromatogram A')
subplot(2,1,2);hist(tB/60); %conversion to minutes
xlabel('retention time (min)')
ylabel('Signal (arbitary units)')
title('chromatogram B')