Code covered by the BSD License  

Highlights from
Lognormal Perfusion Model

image thumbnail

Lognormal Perfusion Model

by

 

13 Mar 2013 (Updated )

Fit disruption replenishment time-intensity data with the lognormal perfusion model

logperf_model(A,mu,sigma,HW,Delay,t,Fast)
%
% The Lognormal Perfusion Model
%
% Reference as:
% Hudson et al. Quantification of flow using ultrasound and microbubbles: 
% a disruption replenishment model based on physical principles. 
% Ultrasound in Medicine & Biology (2009) vol. 35 (12) pp. 2007-20
% 
% or
% 
% Hudson et al. The Lognormal Perfusion Model for Disruption Replenishment 
% Measurements of Blood Flow: In Vivo Validation. 
% Ultrasound in Medicine & Biology (2011) vol. 37 (10) pp. 1571-1578
%
% JM Hudson
% March 2013


% Input parameters:
% A = Vascular cross sectional area
% mu = lognormal mean of the flow speed distribution
% sigma = standard deviation of the flow speed distribution
% Delay = 
% mu 	= mean of lognormal distribution
% sigma	= standard deviation of lognormal distribution
% Delay = Delay due to the separation between the disruption and detection
% beams
% t 	= Time
% Fast  = 
% *****************************************************************
%


function [Signal] = logperf_model(A,mu,sigma,HW,Delay,t,Fast)
%  A = 1
%  mu = 1,
%  sigma =0.5
%  HW=8
%  Delay = 0
%  t=0:0.5:25;
%  Fast =1


    tic
    warning off;
    Amp =1;

    % This determines the resolution of 'D' and determines the speed of the processing
    % Set Fast to 1 for quick fitting and 0 for better plots

    if Fast == 1
        Res = .3;
    else
        Res = 0.15;
    end 

    A = abs(A);
    mu=(mu); 
    HW=abs(HW); 
    Delay=(Delay);

   if HW < 4 || Delay < -1
      D = 0:Res:(2*(10)+Delay);  % Dummy x-axis for simulation     
   else
      D = 0:Res:(2*(HW)+Delay);  % Dummy x-axis for simulation
   end
   
    % Make a big 2D array to calc things quickly
    B = Wsinc(D,HW,Delay,Amp); % Sinc Beam function
    D_big = repmat(D,[length(t) 1]);
    t_big = repmat(t',[1 length(D)]);
    F = A .* 1./2 .* erfc( (log(D_big./t_big) - (mu))./sigma./sqrt(2)); 
    z = F .* repmat(B,[length(t),1]);
    Signal = trapz(D',z');
    if isnan(Signal(1)) == 1
        Signal(1) = 0;
    else
    Signal = Signal-Signal(1);
    end
    

return;



Contact us