Particle tracking using IDL -- John C. Crocker and Eric R. Weeks (this page written primarily by Brenton Hoffman)
Home | Download software | Tutorial | Extra software

Click here to return to microrheology summary page


"Two-point microrheology of inhomogeneous soft materials"
JC Crocker, MT Valentine, ER Weeks, T Gisler, PD Kaplan, AG Yodh, and DA Weitz, Phys. Rev. Lett. 85, 888 (2000).
  • View abstract / PDF version / Web page related to this work

    A program designed to use multiple particle trajectories to isolate the correlated motion between particles (i.e. the motion of one particle due to movement of another). The outputs are the parallel (denoted as rr) and perpendicular (denoted as tt or theta-theta) components of the correlation tensor as a function of lagtime and particle separation distance. These correspond to the diagonal components of a correlation tensor that would be called 'S' in the Fluctuation/Dissipation Theorem in Chaikin & Lubensky. With some further analysis, using, the output can generate shear moduli with

    One advantage of using correlated motion instead of single particle motion is that mechanical properties can be accurately determined even if small scale heterogeneities exist in the material of interest. However, the material being probed must act as a continuum on the length scales probed. The correlation tensor should decay as 1/r in a continuum. Therefore, the spatial behavior of the correlated motion must be probed before moduli data can be accurately determined.


    The required inputs to this program are a trajectory file, either a filename or an idl variable, and parameters to control the program. There are four required parameters. The first is the minimum distance over which to correlate particle motion. At small separation distances particles' images tend to overlap, leading to particle tracking errors. The second parameter sets the maximum distance over which to correlate particle motions. In thin samples, this parameter should be limited to a regime with continuum like behavior. The third parameter is the number of bins the spatial data will be stored in. More bins give more points to evaluate the spatial behavior, but also decrease the statistical power in each point. The choice of this parameter should be optimized accordingly. The fourth parameter is the maximum lagtime over which particles will be correlated. This time should be set to less than the maximum time in the trajectories.


    Note that the first two (MICPERPIX, TIMESTEP) are required for the rheology measurements to have meaningful units.


    IDL> m2=msd2pnt(t,[2,8,5,10],mic=0.096, tim=1./1000,out="msd2.test")
    where t is a trajectory file, 2 and 8 and are the minimum and maximum particle separation distances to be used in microns, 5 is the number of bins in the spatial data, 10 is the maximum time in seconds, the pixel scale is 0.096 micron / pixel, the frame rate is 1000 fps, and the results will be written to msd2.test.


    The result 'm2' has the form:
    	In 2 dimensions (polar coordinates):
     	  res(*,*,0) contains the log-spaced 'dr' values.
    	  res(*,*,1) contains the log-spaced 'dt' values.
    	  res(*,*,2) is the longitudinal, 'r-r' mean component.
    	  res(*,*,3) is the transverse, 'theta-theta' mean component.
    	  res(*,*,4:5) stores the *variances* of the r,theta parts.
    	  res(*,*,6) contains the total number of points used in each bin.
    	In 3 dimensions (spherical coordinates):
    	  res(*,*,0) contains the log-spaced 'dr' values.
    	  res(*,*,1) contains the log-spaced 'dt' values.
    	  res(*,*,2) is the longitudinal, 'r-r' mean component.
    	  res(*,*,3) is the transverse, 'theta-theta' mean component.
    	  res(*,*,4) is the transverse, 'phi-phi' mean component.
    	  res(*,*,5:7) store the *variances* of the r,theta,phi parts.
    	  res(*,*,8) contains the total number of points used in each bin.
    The means are what you want, the variances and numbers are intended for use in error estimation of the mean values. The units are um^2 if the inputs are in um. The spatial variation of the two point correlation tensor can be checked by plotting the longitudinal or transverse components as a function of r, the particle separation. This operation should be done at several time points to ensure continuum like behavior at all lagtimes. An example is shown below:

    IDL> t=10
    IDL> Drr=(m2(*,t,2))
    IDL> r=m2(*,t,0)
    IDL> plot_oo,r,Drr,psym=circ()
    IDL> print,m2(0,t,1); prints the lagtime value this data corresponds to

    This would plot Drr as function of r at the 10th time point. The resulting plot should have a slope of negative one. Different lagtimes can be evaluated by changing the value of t. Different slopes indicate non-continuum behavior. The range of particle separations can be set to something large here and then limited in


    The array structure may appear confusing at first glance -- why is this a 3D array, for example? Think of the output as a short stack of planes. In two dimensions, there are 6 planes. The first plane [m2(*,*,0)] is the values of dr, and is repetitive:
    r1 r1 r1 r1 r1 ...
    r2 r2 r2 r2 r2 ...
    r3 r3 r3 r3 r3 ...
    The second plane [m2(*,*,1)] contains the values of dt, and is also repetitive; but runs at "right angles" to the dr information:
    t1 t2 t3 t4 t5 ...
    t1 t2 t3 t4 t5 ...
    t1 t2 t3 t4 t5 ...

    Thus if you look at the data along the rows in each plane, you are taking the data at constant r, but varying t. Likewise the data along the columns is constant t, but varying r. The correlation functions you want are stored in the higher planes, m2(*,*,2) etc. Thus, you can generate plots as a function of dr or dt:

    IDL> n = 15
    IDL> plot_oo,m2(*,n,0),m2(*,n,2),psym=circ(),xtitle="dr"
    IDL> print,"lag time dt = ",m2(0,n,1)

    IDL> n = 20
    IDL> plot_oo,m2(n,*,1),m2(n,*,2),psym=circ(),xtitle="dt"
    IDL> print,"separation dr = ",m2(n,0,0)

    In the first example, you're plotting the data from plane 2 (the correlation) as a function of plane 0 (the dr values, taken at constant dt). You take a quick peek at plane 1 to see what the dt is that you've chosen by setting n=15. In the second example, you plot the data from plane 2 as a function of plane 1 (the dt values, at constant dr). You then look at plane 0 to see what the dr is that corresponds to n=20.

    Or, you can take the whole thing and feed it into IDL's surface plotter:

    IDL> surface,m2(*,*,2),m2(*,*,0),m2(*,*,1),/xl,/yl,yr=[0.1,10]

    And that is the benefit of having this in a 3D array, it's exactly what IDL's surface command wants to see. But, I personally (Eric) prefer to examine several individual curves rather than the whole surface. For example, it's nice to check if you have 1/r behavior, so I like to plot a few curves (at different dt's) together on a log-log plot to see if they are 1/r.

    This program was written by John, and this webpage was written by Brenton Hoffman, a student working with John (with some editing to html-ize it by Eric). Eric wrote the subsection entitled "Array Structure".

    Click here to return to microrheology summary page

    Contact us