fit_av.c File Reference

Fit A-V difference using compartmental model. More...

Go to the source code of this file.

Functions
int	simC3vb (double *t, double *cab, const int nr, double f, double k1, double k2, double k3, double k4, double k5, double k6, double *cvb, double *ct)

Detailed Description

Fit A-V difference using compartmental model.

Copyright

(c) Turku PET Centre

Author

Vesa Oikonen

Definition in file fit_av.c.

Function Documentation

simC3vb()

int simC3vb	(	double *	t,
		double *	cab,
		const int	nr,
		double	f,
		double	k1,
		double	k2,
		double	k3,
		double	k4,
		double	k5,
		double	k6,
		double *	cvb,
		double *	ct )

Simulate venous blood TAC using 1-3 tissue compartment model (compartments in series)

Memory for cvb must be allocated in the calling program. To retrieve the tissue TAC, pointer to allocated memory for ct can be given.

The units of rate constants must be related to the time unit; 1/min and min, or 1/sec and sec.

Returns

Function returns 0 when successful, else a value >= 1.

Author

Vesa Oikonen

See also

simC3vs, simC1, simC3p, simC3vp

Parameters

t	Array of time values; must be in increasing order.
cab	Array of arterial blood activities.
nr	Number of values (samples) in TACs.
f	Perfusion; must be f>=K1.
k1	Rate constant of the model; must be K1<=f.
k2	Rate constant of the model.
k3	Rate constant of the model.
k4	Rate constant of the model.
k5	Rate constant of the model.
k6	Rate constant of the model.
cvb	Pointer for venous blood TAC array to be simulated; must be allocated
ct	Pointer for tissue TAC to be simulated, or NULL

Definition at line 109 of file fit_av.c.

  {
  int i;
  double b, c, d, w, z, dt2;
  double cai, ca_last, t_last;
  double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
  double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
 
 
  /* Check for data */
  if(nr<2) return 1;
  if(cvb==NULL) return 2;
 
  /* Check actual parameter number */
  if(!(f>=0.0) || !(k1>=0.0) || k1>f) return 3;
  if(k3<=0.0) {k3=0.0; if(k2<=0.0) k2=0.0;}
  else if(k5<=0.0) {k5=0.0; if(k4<=0.0) k4=0.0;}
  else {if(k6<=0.0) k6=0.0;}
 
  /* Calculate curves */
  t_last=0.0; if(t[0]<t_last) t_last=t[0]; 
  cai=ca_last=0.0;
  ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=0.0;
  ct1=ct2=ct3=ct1i=ct2i=ct3i=0.0;
  for(i=0; i<nr; i++) {
    /* delta time / 2 */
    dt2=0.5*(t[i]-t_last);
    /* calculate values */
    if(dt2<0.0) {
      return 5;
    } else if(dt2>0.0) {
      /* arterial integral */
      cai+=(cab[i]+ca_last)*dt2;
      /* partial results */
      b=ct1i_last+dt2*ct1_last;
      c=ct2i_last+dt2*ct2_last;
      d=ct3i_last+dt2*ct3_last;
      w=k4 + k5 - (k5*k6*dt2)/(1.0+k6*dt2);
      z=1.0+w*dt2;
      /* 1st tissue compartment and its integral */
      ct1 = (
          + k1*z*cai + (k3*k4*dt2 - (k2+k3)*z)*b
          + k4*c + k4*k6*dt2*d/(1.0+k6*dt2)
        ) / ( z*(1.0 + dt2*(k2+k3)) - k3*k4*dt2*dt2 );
      ct1i = ct1i_last + dt2*(ct1_last+ct1);
      /* 2nd tissue compartment and its integral */
      ct2 = (k3*ct1i - w*c + k6*d/(1.0+k6*dt2)) / z;
      ct2i = ct2i_last + dt2*(ct2_last+ct2);
      /* 3rd tissue compartment and its integral */
      ct3 = (k5*ct2i - k6*d) / (1.0 + k6*dt2);
      ct3i = ct3i_last + dt2*(ct3_last+ct3);
    }
    if(f>0.0) {
      double dct = k1*cab[i] - k2*ct1; 
      cvb[i] = cab[i] - dct/f;
    } else
      cvb[i]=cab[i];
 
    /* copy values to argument arrays; set very small values to zero */
    if(ct!=NULL) {ct[i]=ct1+ct2+ct3; if(fabs(ct[i])<1.0e-12) ct[i]=0.0;}
 
    /* prepare to the next loop */
    t_last=t[i]; ca_last=cab[i];
    ct1_last=ct1; ct1i_last=ct1i;
    ct2_last=ct2; ct2i_last=ct2i;
    ct3_last=ct3; ct3i_last=ct3i;
  }
 
  return 0;
}