54 lines
2.7 KiB
Matlab
54 lines
2.7 KiB
Matlab
function [LoadingRate, StandardError, ConfidenceInterval] = calculateLoadingRate(this, ovenObj)
|
|
n = this.NumberOfAtoms;
|
|
DynamicalQuantities = this.ParticleDynamicalQuantities;
|
|
CollisionEvents = zeros(1, n);
|
|
NumberOfLoadedAtoms = zeros(1, n);
|
|
|
|
|
|
% Include the stochastic process of background collisions
|
|
for AtomIndex = 1:n
|
|
this.TimeSpentInInteractionRegion(AtomIndex) = this.computeTimeSpentInInteractionRegion(squeeze(DynamicalQuantities(AtomIndex,:,1:3)));
|
|
CollisionEvents(AtomIndex) = this.computeCollisionProbability(ovenObj, this.TimeSpentInInteractionRegion(AtomIndex));
|
|
end
|
|
|
|
% Count the number of loaded atoms subject to conditions
|
|
|
|
switch this.ErrorEstimationMethod
|
|
case 'bootstrap'
|
|
NumberOfTimeSteps = int64(this.SimulationTime/this.TimeStep);
|
|
NumberOfLoadedAtoms = zeros(1, NumberOfTimeSteps);
|
|
LoadedAtomIndices = [];
|
|
for TimeIndex = 1:NumberOfTimeSteps
|
|
if TimeIndex ~= 1
|
|
NumberOfLoadedAtoms(TimeIndex) = NumberOfLoadedAtoms(TimeIndex-1);
|
|
end
|
|
for AtomIndex = 1:n
|
|
Position = squeeze(DynamicalQuantities(AtomIndex, TimeIndex, 1:3))';
|
|
if this.exitCondition(Position, CollisionEvents(AtomIndex))
|
|
if ~ismember(AtomIndex, LoadedAtomIndices)
|
|
NumberOfLoadedAtoms(TimeIndex) = NumberOfLoadedAtoms(TimeIndex) + 1;
|
|
LoadedAtomIndices(end+1) = AtomIndex;
|
|
end
|
|
else
|
|
if ismember(AtomIndex, LoadedAtomIndices)
|
|
NumberOfLoadedAtoms(TimeIndex) = NumberOfLoadedAtoms(TimeIndex) - 1;
|
|
LoadedAtomIndices(LoadedAtomIndices==AtomIndex) = [];
|
|
end
|
|
end
|
|
end
|
|
end
|
|
[LoadingRate, StandardError, ConfidenceInterval] = this.bootstrapErrorEstimation(ovenObj, NumberOfLoadedAtoms);
|
|
case 'jackknife'
|
|
for AtomIndex = 1:n
|
|
if AtomIndex ~= 1
|
|
NumberOfLoadedAtoms(AtomIndex) = NumberOfLoadedAtoms(AtomIndex-1);
|
|
end
|
|
Position = squeeze(DynamicalQuantities(AtomIndex, end, 1:3))';
|
|
if this.exitCondition(Position, CollisionEvents(AtomIndex))
|
|
NumberOfLoadedAtoms(AtomIndex) = NumberOfLoadedAtoms(AtomIndex) + 1;
|
|
end
|
|
end
|
|
[LoadingRate, StandardError, ConfidenceInterval] = jackknifeErrorEstimation(this, ovenObj, NumberOfLoadedAtoms);
|
|
end
|
|
end
|