function ret = calculateLocalSaturationIntensity(~, PeakIntensity, PositionVector, WaveVectorOrigin, WaveVectorEndPoint, BeamRadius, BeamWaist)
    WaveVector                          = WaveVectorEndPoint - WaveVectorOrigin; % Line
    PositionVectorFromWaveVectorOrigin  = PositionVector     - WaveVectorOrigin; % Point = PositionVector
    
    %Height of parallelogram (Distance between point and line) = Area of parallelogram / Base
    %One side of parallelogram = PositionVectorFromWaveVectorOrigin
    %Base = Wavevector
    %Area = One side of parallelogram X Base
    DistanceBetweenAtomAndLaserBeamAxis = norm(cross(PositionVectorFromWaveVectorOrigin, WaveVector))./ norm(WaveVector);
    
    if DistanceBetweenAtomAndLaserBeamAxis <= BeamRadius
        ret = PeakIntensity * exp(-2*DistanceBetweenAtomAndLaserBeamAxis^2 / BeamWaist^2);
    else
        ret = 0;
    end
end