Concept
DigImage: Particle Tracking
:::: IMAGE CAPTURE
:::: PARTICLE LOCATION
:::: PARTICLE MATCHING
:::: OPERATION
Subsequent Analysis
:::: STRUCTURE OF FILES
:::: BASIC OPERATION OF TRK2DVEL
:::: CONFIGURATION OF TRK2DVEL
:::: FACILITIES
Advanced Features
:::: PARTICLE-TAGGED DATA
:::: IMAGE PREPROCESSING
Experimental Hints
Index
Once all the particles in an image have been found (at t = tn+1, say), they need to be related back to the previous image (t = tn, say) to determine which particle image is which physical particle. In DigImage we use a modification of what is known in operations research as the Transportation Algorithm. While the problem solved by the transportation algorithm may be represented as a 0-1 totally unimodular integer linear program, it is more efficient and illuminating to take a graph theory approach.
The idea is to choose a set of associations between two sets of entities, such that the set of associations is optimal in the sense that it minimises some linear function of the associations it includes. For the particle tracking, one of the sets is the set of particles P at t = tn and the other the set of particles Q at t = tn+1. We shall start by assigning a label to all the particles images in the two images. At t = tn the particle images are labelled pi for i=1 to i=M, while at t = tn+1 they are labelled qj for j=1 to j=N. We now define a set of association variables aij. If aij is equal to one, then we will say that pi at t = tn is produced by the same particle as qj at t = tn+1. If aij is zero, then pi and qj represent different physical particles.
For the time being we shall assume that there is one and only one physical particle for each of the particle images. We shall consider groups of particles later in this discussion. For the present it is obvious that, for given pi, at most only one value of j can give aij equal to one, otherwise the physical particle must be two places at once! Identical arguments apply for each pj. If M is equal to N, it may be possible for there to be exactly M = N values of aij equal to one. However, this will seldom happen in real experiments, where there will normally be fewer than M = N values of aij equal to one. Moreover, the number of particles images at the two times will not always be equal.
There are many reasons why the number of particles in the image may be different at t = tn and t = tn+1. The simplest is that the particle may have moved outside the region of the flow being tracked, either by moving outside the bounds of the tracking window, or by moving out of the illuminated region (e.g. moving out of a sheet of light). To overcome this problem we define a0j and ai0 as dummy particles at times t = tn and t = tn+1. Unlike ordinary particles, more than one value of j or i may give a non zero value of a0j or ai0 (respectively). In this case a non zero value of ai0 indicates that particle pi at t = tn has been lost from the image by t = tn+1, either by moving out of the image or for some other reason. Similarly, a0j = 1 represents a particle qj present at t = tn+1 which was not there at t = tn.
In order to determine the optimal set of non zero aij, we must first define the functional to be optimised. The only restriction this method puts on the functional is that it is linear in the associations, aij, and so may be represented by Z, the sum over i and j of aijcij. Elements of cij represent the cost of associating particle pi at t = tn with particle qj at t = tn+1. The optimal solution will be chosen to minimise the objective function Z.
Typically the costs cij will be specified using some function of the particle positions, particle characteristics, temporal history and the physics of the flow. Conceptually the simplest model is to set cij equal to the separation between particle pi and particle qj (c0j and ci0 may be set to the distance to the boundaries of the observed region, or the maximum allowable distance a particle may be allowed to travel between tn and tn+1). The optimal solution will then try to minimise the particle displacements, allowing only associations which do not exceed the cost limits placed by c0j and ci0. The costs cij could equally as easily be the squares of the displacements, yielding a type of least squares optimal solution.
If we are trying to measure the fluid velocity (rather than Brownian motion, say), then a more appropriate set of cost functions would include some fluid dynamics. This may be achieved at the most basic level by predicting the positions the particles at t = tn will have at t = tn+1 using their velocity (and possibly acceleration) at t = tn. The costs cij may then be some function of the separation between the predicted position of pi and the position of qj. If a particle at t = tn has only just entered the image, then we are unlikely to have more than a rough estimate for its velocity and so are unable to predict accurately where it might be at t = tn+1. To enable matchings to still occur to such particles, we must reduce the costs of associations with them and allow matchings over larger distances than for particles for which we have a velocity history (we may also, however, add some fixed cost for this new member). While this may produce some mismatching, the requirement for a much more exact match would not then be satisfied at t = tn+2, and so the mismatch would not continue. During subsequent analysis, if we accept only paths which passed through three or more samples during the tracking phase, then we will eliminate any mismatches due to the less stringent matching requirement for a particle with no velocity history.
Additional factors such as the particle size, intensity, shape or even colour may easily be brought into the costing function. Every added component in a well-chosen functional will increase the probability of a correct matching, but at the expense of increased computation. Fortunately, provided the particle seeding density is not too dense, the extra criteria are unlikely to add significantly to the quality of the results. Experience has shown that the tracking results are relatively insensitive to the exact function used for the costs cij. Any mismatches which arise due to a short coming in the costing procedure will be short lived (they will fail to match on the next step) and may be trapped during the subsequent analysis phase through acceleration checks, for example. DigImage includes information about particle geometry and intensity in a simple multiplicative manner.
For example, suppose the basic cost of an association between particle images pi and qj is simply the separation between the predicted position of pi and the position of qj, viz.
| (1) |
where xi, xj are the positions of pi at tn and qj at tn+1 (respectively), ui is a measure of the velocity of pi at tn, and dt is tn+1tn. The cost of the association, cij, is then related to Bij by
| (2) |
where hi is the new particle discount function, ej the ellipticity premium function, j the size premium function, tj the threshold premium function, Iij the background intensity premium function, and Ni the new particle joining fee. We may use additional information about the neighbourhood of the particles pi and qj by introducing the cross correlation function Xij between subregions of the image centred on pi and qj. The costing function Yi then determines the contribution this provides to the total cost.
The various weighting functions are defined as follows:
| Modifier | Description |
| hi | New particle discount. If pi at tn does not have any velocity history (i.e. it was not matched with any particle at tn1), then it will be eligible for a discount of Maximum_matching_distance/(New_particle_v_error*dt), where the numerator (set through [;USPM Maximum matching distance]) is the maximum distance over which a normal particle association may occur, and New_particle_v_error ([;USM Max velocity error for new paths]) is the maximum velocity error for a new particle. If pi at tn has a velocity history, then Ni is unity and there is no new particle discount. |
| ej | Ellipticity premium. If the ellipticity of pj falls within specified bounds (set by [;USLE Ellipticity limit]), then ej takes the value of one. If qj is outside the bounds, then ej will be infinite if qj is to be ignored. Alternatively, qj may be considered as two particles, the first of which has ej = 1, while the second incurs the premium specified (set by [;USPE Premium if ellipticity exceeded]). |
| j | Size premium. If the size of the particle lies between specified limits (set by [;USLR Lower size limit] and [;USLS Upper size limit]), then j is unity. If qj is too small, then j is either an inflated value (specified by [;USPR Premium if particle too small]) or infinity, depending on whether the particle is to be included or not. If qj is too large, then it may either be ignored (infinite j), or treated as two particles, the first of which has j = 1 and second of which incurs the premium specified (by [;USPS Premium if particle too large]). |
| tj | Threshold premium. As noted earlier, DigImage employs two or more separate thresholds during the particle location phase. Particles qj which satisfy threshold above the midpoint have tj = 1, while those which satisfy only the less stringent (lower) thresholds have j equal to the value specified (in [;USPI Premium if faint particle]). |
| Iij | Background intensity premium. If background intensity pricing is enabled (through [;USPB Enable background intensity pricing]) and 511 or fewer particles are being tracked, then if the local background intensity in the neighbourhood of pi and pj differ by two levels (the range is divided into eight levels as specified by [;USLB: Background intensity evaluation]), the association suffers a price premium (specified by [;USPC Premium if two intensity steps]). If the intensity difference is more than two levels, then the premium is greater (specified by [;USPD Premium if many intensity steps]). On the other hand, if the intensity difference is one or less, or background intensity pricing is not enabled, then Iij = 1. |
| Ni | New particle joining fee. If pi at tn does not have any velocity history (i.e. it was not matched with any particle at tn1), then in addition to the cost based on the separation between its predicted position and the location where qj is found, it will also suffer a new member fee equal to Ni. Typically the joining fee will be approximately half the maximum cost that may be imposed on the association. If pi does have a velocity history, then there is no joining fee and Ni is zero. The new particle joining fee is set by [;USPJ New particle joining fee]. In most cases the default value of 0.25 is appropriate. |
| Yi | Cross correlation costing. This factor multiplies 1 Xij to determine how the cross correlation information is costed in the solution process. Two user-specified costings are implemented, one for particles with no velocity history (set using [;USPK New particle correlation costing]) and a separate one for particles with a velocity history (set using [;USPL Particle correlation costing with history]). Note that correlation costing is relatively expensive to compute. In most situations it adds little to the accuracy of the tracking process, and so may be disabled by setting Yi to zero. This is the default value. |
Note that this tracking technique is able to cope with flows which have velocity gradients parallel to the viewing direction. As spatial correlation is not used, there is no confusion if a particle at one side of a light sheet (say) is travelling in the opposite direction to one on the other side, even if their paths appear to cross on the video image.
For particles with no velocity history, DigImage will make an
estimate of their velocity prior to evaluating Bij.
The [;USN: New particle behaviour]
menu controls how this estimate is made:
| Access sequence | Description |
| [;USNC Constant predicted velocity] | Assumes on average the particles are moving with some spatially uniform velocity. |
| [;USNL Local velocity estimate] | Determines the local mean velocity from the other particles and uses this to estimate the velocity for the new particle. |
| [;USNR Region-dependent predicted velocity] | Assumes on average the particles within one region are moving with some spatially uniform velocity, and those in the second region are moving with a different spatially uniform velocity. |
| [;USNZ Zero predicted velocity] | Assumes on average the particles have a zero mean velocity. |
The method by which the optimal set of associations is found closely resembles the transportation algorithm. Here we outline the principle of the solution algorithm rather than the specific implementation used by DigImage.
Once the costs cij have been specified, any value of cij exceeding c0j or ci0 (i.e. the cost of a particle leaving the region which is given by the lesser of the cost based on its distance to the boundary and that based on the maximum matching distance [;USPM Maximum matching distance]) and the cost is set to infinity: that matching will never occur so we need not consider it. An initial guess is made at aij such that, for given i = I, j = J is chosen such that cIJ is the minimum of cIj, and aIJ set to one, provided no value of aiJ is already one. Frequently this initial solution will be very close to optimal, but that will not always be the case. It is now necessary to iterate until the optimal solution is found.
For each iteration, we scan through the list of zero aij which have finite costs, and evaluate the cost of bringing that association into the solution. Consider the zero association aIJ and suppose that aIJ = 1 and aKJ = 1. Now if aIJ were to be brought in, then aIL and aKJ must leave the solution, otherwise I and J will be two places at once. Further, as particle images i = K and j = L must be related to some physical particle, it will be necessary to let aKL enter the solution. Thus the reduced cost of bringing in aIJ is cIJ+cKLcILcKJ. If this reduced cost is negative, then bringing aIJ into the solution is favourable (will decrease the objective function). After scanning through all zero aij with finite costs, we bring the aij in which had the most negative reduced cost and start the next iteration. If all aij have reduced costs greater than or equal to zero, then the optimal solution has been found.
This matching technique is executed twice for every time step in DigImage. The first time, we match t = tn with t = tn+1. Any particles pi at t = tn for which we do not find a match at t = tn+1 is stored for later use. Particles qj at t = tn+1 for which there was no match at t = tn are then compared with particles at t = tn1 for which there was no match at t = tn (provided there was a velocity history at t = tn1). Thus we are able to follow particle paths even if the particle concerned is not visible for one sample. During the subsequent analysis stage care should be exercised in using any data which spans such a disappearance; they may be excluded from such analyses. The main reason behind this second matching phase is to maintain a velocity history for that particle so that its future position may be predicted.
Difficulties with particle paths crossing
such that two particles are in the same place at t = tn+1
are overcome in a simple manner. We assume that any such conglomeration
of particles will have unusual shape and or size characteristics.
The combined particles may appear like a very large, or very elliptical
(if the two particles are side-by-side) particle. If such an anomalous
particle is detected, then DigImage will treat it as two
particles. The first of the pair is treated in the ordinary
manner, with the association costs prescribed in the manner outlined
above. The second will be treated with more caution as
it may not be a second particle. The cost of associating it with
any particle at t = tn
will be increased by some predefined factor ([;USPE
Premium if ellipticity exceeded]
or [;USPR Premium if particle too
large]) depending on why
the particle is to be considered as two. These premiums are related
to the probability that the anomalous blob is infect two particles.
At present these factors represent a linear increase in the cost
as was outlined above (ejj).
In practice these price premiums are seldom invoked during the
particle tracking as the relatively sparse particle seeding ensures
the probability of two particles coinciding is low. Moreover,
if they are invoked incorrectly, then it is unlikely that the
path containing the erroneous matching will survive. Note that
matching two samples back (i.e. to tn1)
is only permitted on the normal part of such a particle.
Goto next document (Operation)
DigImage documentation page
DigImage home page
Stuart Dalziel's home page