In the accompanying drawings:
a-d illustrates an example of the estimation of target position, lateral velocity, and road curvature parameters for a straight roadway;
a-b illustrate an example of the target state RMS errors from unconstrained and constrained filtering on the straight roadway, corresponding to
a-d illustrate an example of the estimation of target position, lateral velocity, and road curvature parameters for a curved roadway;
a-b illustrate an example of the target state RMS errors from unconstrained and constrained filtering for the curved roadway, corresponding to
a-d illustrate an example of the estimation of target position, lateral velocity, and associated RMS errors for a straight roadway involving a lane change;
a-d illustrates an example of the estimation of target position, lateral velocity, and their RMS errors for a curved roadway involving a lane change;
a illustrates a geometry of a bicycle model of a vehicle undergoing a turn;
b illustrates a geometry of the steered wheel illustrated in
Referring to
Referring to
Referring to
Referring also to
A well-designed and constructed roadway 34 can be described by a set of parameters, including curvature, wherein the curvature of a segment of the roadway 34 is defined as:
where R is the radius of the segment. In general, for a piece of smooth roadway 34, the curvature variation can be described as a function of a distance l along the roadway 34 by a so-called clothoid model, i.e.:
where C1=1/A2 and A is referred to as the clothoid parameter.
Referring to
Substituting equation (2) into equation (3) gives
Δθ=θ−θ0=C0l+C1l2/2 (4)
Referring to
Assuming the heading angle θ to be within 15 degrees, i.e. |θ|<15°, equations (5) and (6) can be approximated by:
Δx=x−x0≈l (7)
Accordingly, the roadway 34 is modeled by an incremental road equation in terms of curvature coefficients (or parameters): C0 and C1. This incremental road equation describes a broad range of road shapes as follows: 1) Straight roadway 34: C0=0 and C1=0; 2) circular roadway 34: C1=0; and 3) a general roadway 34 with an arbitrary shape for which the change in heading angle θ is less than 15 degrees: C0>0.
The road curvature parameters C0 and C1 are estimated using data from motion sensors (yaw rate sensor 16 and speed sensor 18) in the host vehicle 12, based upon the assumption that the host vehicle 12 moves along the center line 41 of the roadway 34 or associated host lane 38.
The road curvature parameters C0 and C1 can be calculated from data of ω, {dot over (ω)}, U, {dot over (U)} responsive to measurements of yaw rate ωh and speed Uh of the host vehicle 12 from the available host vehicle 12 motion sensors. However, generally the measurements of yaw rate ωh and speed Uh, from the yaw rate sensor 16 and speed sensor 18 respectively, are noisy. A host state filter implemented by a first Kalman filter 52 is beneficial to generate estimates of ω, {dot over (ω)}, U, {dot over (U)} from the associated noisy measurements of yaw rate ωh and speed Uh; after which a curvature filter implemented by a second Kalman filter 54 is used to generate smoothed estimates of the curvature parameters C0 and C1. The dynamics of the host vehicle 12 for the host state filter follows a predefined set of kinematic equations (constant velocity in this case) given by:
xk+1h=Fkh·xkh+wkh, wkh˜N(0, Qkh) (9)
zkh=Hkh·xkh+vkh, vkh˜N(0, Rkh) (10)
where
and where T is the sampling period, superscript (•)h is used to indicate that the filter is the host filter, and Uh and ωh are host vehicle 12 speed and yaw rate measurements. The first Kalman filter 52 is implemented to estimate the host state {circumflex over (x)}k|kh and its error covariance as illustrated in
The estimate of the host state from the first Kalman filter 52, i.e. the host state filter, is then used to generate a synthetic measurement that is input to the second Kalman filter 54, i.e. curvature coefficient (or parameter) filter, wherein the associated Kalman filters 52, 54 operate in accordance with the Kalman filtering process described more fully in the Appendix hereinbelow. The relationship between the road curvature parameters C0, C1 and the host state variables ω, {dot over (ω)}, U, {dot over (U)} is derived as follows:
From equation (4), the radius R of road curvature is expressed generally as a function R(l) of the distance l along the roadway, as is illustrated in
{dot over (θ)}=C0·{dot over (l)}+C1·l·{dot over (l)}=(C0+C1·l)·{dot over (l)}. (12)
Noting that {dot over (θ)}=ω, the yaw rate of the host vehicle 12, and that {dot over (l)}=U, the speed of the host vehicle 12, and substituting the clothoid model of equation (2) in equation (12), yields:
ω=C·U (13)
or
Clothoid parameter C0 is given as the value of curvature C at l=0, or
Taking the derivative on both sides of equation (14) yields
Using the definition of C1, from equation (2), C1 may be expressed in terms of the host state as follows:
The system equations for the second Kalman filter 54, i.e. the curvature filter, that generates curvature estimates Ĉ0
xk+1C=FkC·xkC+wkC, wkC˜N(0, QkC) (18)
zkC=HC·xkC+vkC, vkC˜N(0, RkC) (19)
where
Δt is the update time period of the second Kalman filter 54, and the values of the elements of the measurement vector zkC are given by the corresponding values of the state variables—i.e. the clothoid parameters C0 and C1—of the curvature filter.
The measurement, zkC, is transformed from the estimated state [Û, {dot over (Û)}, {circumflex over (ω)}, {dot over ({circumflex over (ω)})}]kT as follows:
and the associated covariance of the measurements is given by:
RkC=JkCPk|kh(JkC)T (22)
where
It should be understood that other systems and methods for estimating the curvature parameters of the roadway 34 may be substituted in the road curvature estimation subsystem 42 for that described above. For example, the curvature parameters of the roadway may also be estimated from images of the roadway 34 by a vision system, either instead of or in conjunction with the above described system based upon measurements of speed Uh and yaw rate ωh from associated motion sensors. Furthermore, it should be understood that yaw rate can be either measured or determined in a variety of ways, or using a variety of means, for example, but not limited to, using a yaw gyro sensor, a steering angle sensor, a differential wheel speed sensor, or a GPS-based sensor; a combination thereof; or functions of measurements therefrom (e.g. a function of, inter alia, steering angle rate).
Referring again to
The result from the coordinate transformation in step (512) of the output from the extended Kalman filter 56 is then partitioned into the following parts, corresponding respectively to the x and y position of the target vehicle 36 relative to the host vehicle 12, wherein the superscript 1 refers to the unconstrained target state of the target vehicle 36:
Referring again to
Prior to discussing the process of steps (514) through (524) for determining whether the target is likely constrained by a constraint, and if so, what is the most likely constraint, the process of fusing the unconstrained target state with state of a constraint will first be described for the case of a target vehicle 36 moving in the same lane as the host vehicle 12. The constraints are applied in the y-direction and are derived from road equations where y-direction state variables are functions of x-direction state variables, consistent with the assumptions that the host vehicle 12 moves along the center line 41 of its lane 38 steadily without in-lane wandering and that the road curvatures of all the parallel lanes 38, 40 are the same, and given that the absolute coordinate system is fixed on the host vehicle 12 at the current instant of time. Assuming the target vehicle 36 is moving in the same lane 38 as the host vehicle 12, and using the road constraint equation with the estimated coefficients (or parameters), in step (514), the constraint state variables are then given in terms of the lateral kinematic variable as:
In step (528), the two y-coordinate estimates, one from the main filter and the other from the road constraint, are then fused as follows:
Finally, the composed estimate of the target state is
where
In step (530), this composed estimate would then be output as the estimate of the target state if the target vehicle 36 were to be determined from steps (514) through (524) to be traveling in the host lane 38.
Returning to the process of steps (514) through (524) for determining whether the target is likely constrained by a constraint, and if so, what is the most likely constraint; according to the assumption that targets follow the same roadway 34, if the target vehicle 36 were known to travel in a particular lane, it would desirable to use estimated road parameters for that lane as a constraint in the main filter of estimating target kinematics. However, the knowledge of which lane the target vehicle 36 is current in is generally not available, especially when the target is moving on a curved roadway 34. Since the road equation (8) is only for the host lane 38 in the host-centered coordinate system, constrained filtering would require knowing which lane the target is in, and different constraint equations would be needed for different lanes. Ignoring the difference of road curvature parameters among these parallel lanes, i.e. assuming the curvature of each lane to be the same, the road equation for an arbitrary lane can be written as:
where B is the width of the lanes and m represents the lane to be described (m=0 corresponds the host lane 38, m=1 corresponds the right neighboring lane 40, m=−1 corresponds the left neighboring lane 40, and so on). Without the prior knowledge of the target lane position, each of the multiple constraints forming a multiple constraint system (analogous to the so-called multiple model system) is tested to determine which, if any, of the constraints are active. A multiple constraint (MC) system is subjected to one of a finite number NC of constraints. Only one constraint can be in effect at any given time. Such systems are referred to as hybrid—they have both continuous (noise) state variables as well as discrete number of constraints.
The following definitions and modeling assumptions are made to facilitate the solution of this problem:
Constraint equations:
yt
where ft
Constraint: among the possible NC constraints
ft
t
j: state estimate at time tk using constraint ft
yt
j,
μt
Constraint jump process: is a Markov chain with known transition probabilities
P{ft
To implement the Markov model—for systems with more than one possible constraint state—it is assumed that at each scan time there is a probability pij that the target will make the transition from constraint state i to state j. These probabilities are assumed to be known a priori and can be expressed in the probability transition matrix as shown below.
The prior probability that fj is correct (fj is in effect) is
P(fj|Z0)=μt
where Z0 is the prior information and
since the correct constraint is among the assumed NC possible constraints.
The constrained target state estimation subsystem 46 provides for determining whether the target state corresponds to a possible constrained state, and if so, then provides for determining the most likely constrained state.
A multiple constraint (MC) estimation algorithm mixes and updates NC constraint-conditioned state estimates using the unconstrained state estimate ŷt
1. Estimation of state variables from multiple constraints: In step (514), using the multiple lane road equation (34) to replace the first row in equation (25), the multiple constraint state estimates are given by:
where
and B is the width of a lane. Stated in another way, the constraint state estimates corresponds to—e.g. matches—the y locations of the centerlines of each possible lane in which the target vehicle 36 could be located.
The associated covariance is given by:
where Ak1 and Ak2 are given by equation (27) and equation (28), Pxt
2. Constraint-conditioned updating: In step (516), the state estimates and covariance conditioned on a constraint being in effect are updated, as well as the constraint likelihood function, for each of the constraints j=1, . . . NC. The updated state estimate and covariances corresponding to constraint j are obtained using measurement ŷt
3. Likelihood calculation: In step (518), the likelihood function corresponding to constraint s is evaluated at the value yt
Λt
wherein the Gaussian distribution N (; ,) has a mean value of
4. Constraint probability evaluations: In step (520), the updated constraint probabilities are calculated for each of the constraints j=1, . . . NC, as follows:
where āj, the probability after transition that constraint j is in effect, is given by
and the normalizing constant is
5. Overall state estimate and covariance: In step (522), the combination of the latest constraint-conditioned state estimates and covariances is given by:
The output of the estimator from step (522) in the above algorithm is then used as the constrained estimates in the fusion process described by equations (29) and (30), and the result of equation (52), instead of the result of equation (33), is used in equation (32).
When the target vehicle 36 is not following the roadway 34 or is changing lanes, imposing the road constraint on target kinematic state variables will result in incorrect estimates that would be worse than using the associated unconstrained estimates. However, noise related estimation errors might cause a correct road constraint to appear invalid. Accordingly, it is beneficial to incorporate a means that can keep the constraints in effect when they are valid, e.g. when the target vehicle 36 follows a particular lane; and lift them off promptly when they are invalid, e.g. when the target vehicle 36 departs from its lane. The unconstrained target state estimate plays a useful role in road constraint validation, since it provides independent target state estimates.
One approach is to test the hypothesis that the unconstrained target state estimate satisfies the road constraint equation, or equivalently, that the constrained estimate and the unconstrained estimate each correspond to the same target. The optimal test would require using all available target state estimates in history through time tk and is generally not practical. A practical approach is the sequential hypothesis testing in which the test is carried out based on the most recent state estimates only. In accordance with the notation used hereinabove, the difference between the constrained and unconstrained target state estimates (y direction only) is denoted:
{circumflex over (δ)}t
as the estimate of
where yt
H0: δt
vs.
H1: δt
The main filter error
{tilde over (y)}t
is assumed independent of the error
which is from the constraints. The covariance of the difference {circumflex over (δ)}t
Assuming that the estimation errors are Gaussian, the test of H0 vs. H1 is as follows:
Accept H0 if ρt
The threshold is chosen such that
P(ρt
where α is a predefined error tolerance value. Note that based on the above Gaussian error assumption, ρt
Based on the above analysis, the hypothesis testing scheme efficiently uses different threshold values for targets in different lanes, with the multiple constraint filtering algorithm providing the knowledge of which lane the target is most likely in currently. Assuming that there are NC possible lanes on the roadway 34, and each lane is described by a constraint equation, the constraint equation with the highest probability μt
The difference between the unconstrained state estimates and lane lt constrained state estimates (y direction only), denoted as:
{circumflex over (δ)}t
is the estimate of
δt
where yt
The test for the “same target” hypothesis is then given by:
H0: δt
vs.
H1: δt
The constrained estimation error is given by:
Assuming that the estimation errors are independent and Gaussian, the test of H0 vs. H1 becomes:
Accept H0 if ρt
where
and the threshold is such that
P(ρt
where
γl
Such a lane adaptive hypothesis testing scheme provides for a prompt switch of the target state estimation output to the unconstrained estimate when the target vehicle 36 leaves its current lane, while the estimation accuracy of a target in host lane 38 is substantially improved by constrained filtering.
In another embodiment of the multiple constraint (MC) estimation algorithm, the constrained state estimate used for the hypothesis testing is the most likely of the separate constrained target state estimates (i.e. in accordance with a “winner take all” strategy), rather than a composite combination of all of the constrained target state estimates. If this most likely constrained state estimate is valid, i.e. if the most likely constrained state estimate corresponds to—e.g. matches—the unconstrained state estimate, then the target state is given by fusing the most likely constrained state estimate and the unconstrained state estimate; otherwise the target state is given by the unconstrained state estimate.
In yet another embodiment of the multiple constraint (MC) estimation algorithm, hypothesis tests are made for each of the constrained state estimates. If none of the hypotheses are accepted, then the target state is given by the unconstrained state estimate. If one of the hypotheses is accepted, then the target state is given by fusing the corresponding constrained state estimate and the unconstrained state estimate. If more than one hypotheses are accepted, then the most likely constrained state may be identified by voting results from a plurality of approaches, or by repeating the hypothesis tests with different associated thresholds.
Generally, the number of constraints (i.e. the number of roadway lanes) can vary with respect to time, as can associated parameters therewith, for example, the width of the lanes of the roadway, so as to accommodate changes in the environment of the host vehicle 12. For example, the host vehicle 12 in one trip could travel on a one-lane road, a two-lane road with opposing traffic, a three-lane road with a center turn lane, a four line road two lanes of opposing traffic, or on a multi-lane divided freeway.
Road vehicle tracking simulations using constrained and unconstrained filtering were carried out for four scenarios. In all scenarios, the host vehicle 12 was moving at 15.5 m/s and a target vehicle 36 is approaching on the same roadway 34 at a speed of 15.5 m/s. The initial position of the target was 125 meters away from the host in the x direction, and the lane width for all lanes was assumed to be 3.6 meters. The measurement variance of the vehicle speed sensor was 0.02 m/s and the variance of the gyroscope yaw rate measurement was 0.0063 rad/s. The variances of radar range, range rate and azimuth angle measurements were 0.5 m, 1 m/s, and 1.5° respectively. Simulation results were then generated from 100 Monte-Carlo runs of the associated tracking filters.
In the first scenario, the host vehicle 12 and the target vehicle 36 were moving on a straight roadway 34 (C0=0 and C1=0) and the target vehicle 36 was moving toward the host vehicle 12 in the same lane.
In the second scenario, the host vehicle 12 and the target vehicle 36 were moving on a curved roadway 34 (C0=−10−5 and C1=−3×10−5) and the target vehicle 36 was moving toward the host vehicle 12 in the same lane.
In the third scenario, the host vehicle 12 and the target vehicle 36 were moving on a straight roadway 34 (C0=0 and C1=0) and the target vehicle 36 was initially approaching in the left neighboring lane. At t=2.2 second (55 radar scans), the target vehicle 36 began to diverge from its lane and turns toward the host lane 38, which resulted in a collision at t=4 seconds (100 radar scans).
The fourth scenario was similar to the third scenario, the only difference being that the vehicles were on a curved roadway 34 (C0=−10−5 and C1=−3×10−5) instead of a straight one. The target vehicle 36 began to diverge at t=2.2 s and results in a collision at t=4 s.
Accordingly, simulation results of road vehicle tracking on both straight and curved roadways 34 show that the predictive collision sensing system 10 could substantially reduce the estimation errors in target vehicle 36 lateral kinematics when the target vehicles 36 were in the host lane 38. When a target vehicle 36 maneuvers from a neighboring lane into the host lane 38, the predictive collision sensing system 10 promptly detects this maneuver and lifts off the road constraint to avoid an otherwise incorrect constrained result. In view of the fact that poor radar angular resolution often results in poor lateral kinematics estimation, the predictive collision sensing system 10 has provided for a substantial improvement in estimation accuracy of target vehicle 36 lateral kinematics, which is beneficial for an early and reliable road vehicle collision prediction.
Referring to
Fy=Cα·α (72)
wherein the associated proportionality constant, Cα—also known as cornering stiffness—is defined as the slope of cornering force Fy with respect to slip angle α at α=0.
Generally, for a vehicle 12 with two laterally displaced steered wheels 60, the associated steer angles δ of different steered wheels 60—i.e. inside and outside relative to a turn—will be different. However, referring to
For a vehicle 12 at a longitudinal speed U following a curved path having a turn radius R, the sum of the lateral forces, i.e. cornering forces Fy, is equal to the product of the mass M of the vehicle 12 times the resulting centripetal acceleration, as follows:
where Fyf and Fyr are the lateral forces at the front 63 and rear 65 wheels respectively. Assuming that the yaw rotational acceleration of the vehicle 12 about the center of gravity CG is negligible, the sum of the moments caused by the front and rear lateral forces is equal to zero, with the result that
Fyf=Fyr·c/b (74)
where b and c are the distances from the center of gravity CG to the front 63 and rear 65 wheels respectively.
The lateral force Fyr at the rear wheel 65 is then given as follows by substituting equation (74) into equation (73):
where b and c are the distances from the center of gravity CG to the front 63 and rear 65 wheels respectively, and Wr is the weight of the vehicle 12 carried by the rear wheel 65. Accordingly, the lateral force Fyr at the rear wheel 65 is given by the product of the portion of vehicle mass (Wr/g) carried by the rear wheel 65 times the lateral acceleration at the rear wheel 65.
Similarly, the lateral force Fyf at the front wheel 63 is given by the product of the portion of vehicle mass (Wf/g) carried by the front wheel 63 times the lateral acceleration at the front wheel 63 as follows:
where Wf is the weight of the vehicle 12 carried by the front wheel 63.
Given the lateral forces Fyf, Fyr at the front 63 and rear 65 wheels respectively, the associated slip angles αf, αr are given from equations (72), (75) and (76) as follows:
αf=Wf·U2/(Cαf·g·R) (77)
and
αr=Wr·U2/(Cαr·g·R) (78)
From the geometry illustrated in
δ=L/R+αf−αr (79)
Substituting for αf and αr in equation (79) from equations (77) and (78) gives:
which can also be expressed as:
δ=L/R+K·αy (81)
where
δ=Steer angle at the front wheels (rad)
L=Wheelbase (m)
R=Turn Radius (m)
U=Longitudinal speed (m/sec)
g=Gravitational acceleration constant=9.81 m/sec2
Wf=Load on the front axle (kg)
Wr=Load on the rear axle (kg)
Cαf=Cornering stiffness of the front tires (kgy/rad)
Cαr=Cornering stiffness of the rear tires (kgy/rad)
K=understeer gradient (rad/g)
ay=lateral acceleration (g)
Equations (80) and (81) describe the relationship between steer angle δ and lateral acceleration ay=U2/(gR). The factor K=[Wf/Cαf−Wr/Cαr]—which provides the sensitivity of steer angle δ to lateral acceleration ay, and which is also referred to as an understeer gradient—consists of two terms, each of which is the ratio of the load on the wheel Wf, Wr (front or rear) to the corresponding cornering stiffness Cαf, Cαr of the associated tires. Depending upon the value of the understeer gradient K, the cornering behavior of the vehicle 12 is classified as either neutral steer, understeer or oversteer, depending upon whether K is equal to zero, greater than zero, or less than zero, respectively.
For a vehicle 12 exhibiting neutral steer,
Wf/Cαf=Wr/Cαr→K=0→αf=αr (82a)
so that, for a constant-radius turn, there would be no change in steer angle δ as the longitudinal speed U is varied.
For a vehicle 12 exhibiting understeer,
Wf/Cαf>Wr/Cαr→K>0→αf>αr (82b)
so that, for a constant-radius turn, the steer angle δ would need to increase with increasing longitudinal speed U in proportion to the product of the understeer gradient K times the lateral acceleration ay.
For a vehicle 12 exhibiting oversteer,
Wf/Cαf<Wr/Cαr→K<0→αf<αr (82c)
so that, for a constant-radius turn, the steer angle δ would need to decrease with increasing longitudinal speed U in proportion to the product of the understeer gradient K times the lateral acceleration ay.
A vehicle 12 steered with a steer angle δ develops a yaw rate ω that is related to the longitudinal speed U and turn radius R by:
Solving for the turn radius R from equation (79) and substituting in equation (83) gives the following relationship between yaw rate ω and steer angle δ:
which can be used to find the relationship between the associated error variances, e.g. assuming that the vehicle 12 exhibits neutral steer behavior. For example, in neutral steer case, K=0, so that equation (84) becomes:
The error variance of ω is given by:
Assuming that the longitudinal speed U and steer angle δ are independent, then
From equation (85) for a neutral steer condition
which when substituted into equation (87) gives:
where σU and σδ are the error variances of longitudinal speed U and steer angle δ.
For a constant turn radius R, from equation (2), C=C0, and from equations (13) and (88) for a neutral steer condition,
from which the relation between curvature error variance σC02 and steer angle error variance σδ2 is given by:
The steer angle sensor 58 can be implemented in various ways, including, but not limited to, an angular position sensor—e.g. shaft encoder, rotary potentiometer or rotary transformer/syncro—adapted to measure the rotation of the steering wheel shaft or input to a steering box, e.g. a pinion of a rack-and-pinion steering box; or a linear position sensor adapted to measure the position of the rack of the rack-and-pinion steering box. The steer angle sensor 58 could be shared with another vehicle control system, e.g. a road following or suspension control system. The steer angle sensor 58 can be used to supplement a yaw rate sensor 16, and can beneficially provide independent information about vehicle maneuvers. Furthermore, the steer angle δ measurement error is substantially independent of longitudinal speed U, in comparison with a gyroscopic yaw rate sensor 16 for which the associated yaw rate ω measurement error is related to vehicle speed, notwithstanding that a gyroscopic yaw rate sensor 16 is generally more accurate and more sensitive to vehicle maneuvers than a steer angle sensor 58 when each is used to generate a measure of yaw angle.
The curvature error variance associated with steer angle δ measurements can be compared with that associated with yaw rate ω measurements in order to identify conditions under which one measurement is more accurate than the other. The error variance of yaw rate ω measured with a gyroscopic yaw rate sensor 16 is given as follows:
σω2=E[(ωm−bm−ω+b)2]=E[(ωm−ω)2+(bm−b)2+2(ωm−ω)(bm−b)] (92)
σω2=E[(ωm−ω)2]+E[(bm−b)2]+2E[(ωm−ω)(bm−b)] (93)
σω2=σωm2+σb2 (94)
where ω is the true yaw rate, ωm is the yaw rate measurement, b is the gyro bias with drift, and bm is the mean gyro bias.
The curvature error variance σc02 of the yaw rate ω is given by equation (97), described hereinbelow. By equating equations (91) and (97), and substituting for σω2 from equation (94), the curvature error variance associated with the steer angle δ measurement is equal to the curvature error variance associated with the yaw rate ω measurement when:
Equation (95) defines a switching curve—e.g. as illustrated in
The error variances and covariance of the road curvature parameters C0 and C1 used by the associated second Kalman filter 54, i.e. the curvature filter, of the road curvature estimation subsystem 42 generally depend upon the quantity 1/U, where U is the longitudinal speed U of the vehicle 12. If U is a Gaussian random variable, analytic solutions for the exact mean and variance of 1/U are not presently known. Instead, U can be assumed to be a non-random variable because the variance σU2 of U is substantially smaller than U in substantially most cases. Accordingly, the various variance and covariance components of the road curvature parameters C0 and C1 can be derived as follows:
Referring to
Generally, the associated state and variance output of the host state filter 52.1 is processed by a curvature estimator 70 so as to provide estimates of the road curvature parameters C0 and C1 and the error covariance associated therewith. Generally, the curvature estimator 70 comprises a first curvature processor 72 which transforms the associated state and covariance output of the host state filter 52.1 to either road curvature parameters C0 and C1, or another related form—comprising a measurement vector Zkc and an associated covariance matrix Rkc thereof—that is either used directly as the output of the road curvature estimation subsystem 42, or is input to a second Kalman filter 54 of the curvature filter 54.1, the output of which is either used as the output of the road curvature estimation subsystem 42, or which is transformed to road curvature parameters C0 and C1 and the associated covariance thereof using a second curvature processor 74. For example, the first 72 and second 74 curvature processors and the host state filter 52.1 of the curvature estimator 70 can be embodied in the signal processor 26.
In accordance with a first embodiment, the curvature estimator 70 comprises the first curvature processor 72 and the second Kalman filter 54, wherein the first curvature processor 72 calculates the road curvature parameters C0 and C1 from the host state [ω, {dot over (ω)}, U, {dot over (U)}]T for input as a measurement vector ZkC to the second Kalman filter 54, e.g. in accordance with equation (21). Similarly, the first curvature processor 72 calculates the associated covariance RkC of the measurement vector from the covariance Pk|kh of the host state vector, e.g. in accordance with equations (21) and (22). The associated second Kalman filter 54 is illustrated in
wherein q=6×10−7, Tk is the sampling time interval, and the variance and covariance elements of the associated R matrix are given by equations (97), (102) and (109). The output of the road curvature estimation subsystem 42 is then given by the output of the second Kalman filter 54. The first embodiment of the curvature estimator 70 is also illustrated in
A second embodiment of the curvature estimator 70 is a modification of the first embodiment, wherein the second Kalman filter 54 is adapted to incorporate a sliding window in the associated filtering process. The length of the sliding window is adapted so as to avoid excessive delay caused by window processing, and for example, in one embodiment, comprises about 5 samples. The associated vectors and matrices—referenced in the Appendix—of the associated second Kalman filter 54 are given by:
where L is the sliding window length,
and q=2×10−7.
The variances and covariance of the associated R matrix are given by equations (97), (102) and (109).
A third embodiment of the curvature estimator 70 is a modification of the second embodiment, wherein the length L of the sliding window is adaptive. For example, the window length can be adapted to be responsive to the road curvature parameters C0 and C1, for example, in accordance with the following rule:
Lk=min{max{25−floor [(|Ĉ0
This rule provides for a larger window length L—e.g. as large as 25 samples—when both C0 and C1 are relatively small, for example, corresponding to a straight section of road. The window length L becomes smaller—e.g. as small as 1 sample—when either C0 or C1 is large, corresponding to a turn or transition of road. Furthermore, the window length L can be sharply reduced to account for a sudden vehicle maneuver, with a limitation on the subsequent increase in window length L to one sample per step so that the previous samples so not adversely affect the output from the curvature estimator 70 when maneuver ends. As the window length L is changed, the associated F, Q, and R matrices in the second Kalman filter 54 are also changed.
In accordance with a fourth embodiment, the curvature estimator 70 comprises the first curvature processor 72, the second Kalman filter 54, and the second curvature processor 74, wherein the first curvature processor 72 calculates the road curvature parameter C0 and Ċ from the host state [ω, {dot over (ω)}, U, {dot over (U)}]T for input as a measurement vector ZkC to the second Kalman filter 54, and the second curvature processor 74 transforms the output of the second Kalman filter 54 to C0 and C1=Ċ/Û as the curvature estimates of the curvature estimator 70. The associated second Kalman filter 54 is illustrated in
wherein q=0.01, Tk is the sampling time interval, and the variance and covariance elements of the associated R matrix are given by equations (97), (106) and (112).
A fifth embodiment of the curvature estimator 70 is a modification of the fourth embodiment, wherein the second Kalman filter 54 is adapted to incorporate a sliding window in the associated filtering process. The length of the sliding window is adapted so as to avoid excessive delay caused by window processing, and for example, in one embodiment, comprises about 5 samples. The associated vectors and matrices—referenced in the Appendix—of the associated second Kalman filter 54 are given by:
wherein L is the sliding window length,
The variances and covariance of the associated R matrix are given by equations (97), (106) and (112).
A sixth embodiment of the curvature estimator 70 is a modification of the fifth embodiment, wherein the length L of the sliding window is adaptive. For example, the window length can be adapted to be responsive to the road curvature parameters C0 and C1, for example, in accordance with the following rule of equation (117). This rule provides for a larger window length L—e.g. as large as 25 samples—when both C0 and C1 are relatively small, for example, corresponding to a straight section of road. The window length L becomes smaller—e.g. as small as 1 sample—when either C0 or C1 is large, corresponding to a turn or transaction of road. Furthermore, the window length L can be sharply reduced to account for a sudden vehicle maneuver, with a limitation on the subsequent increase in window length L to one sample per step so that the previous samples so not adversely affect the output from the curvature estimator 70 when maneuver ends. As the window length L is changed, the associated F, Q, and R matrices in the second Kalman filter 54 are also changed.
In accordance with a seventh embodiment, the curvature estimator 70 comprises the first curvature processor 72, which calculates the road curvature parameters C0 and C1 from the host state [ω, {dot over (ω)}, U, {dot over (U)}]T as the output of the road curvature estimation subsystem 42—without using a second Kalman filter 54.
In accordance with an eighth embodiment, the curvature estimator 70 comprises the first curvature processor 72, which determines the road curvature parameters C0 and C1 of the clothoid model by a curve fit of a trajectory of the host state [ω, {dot over (ω)}, U, {dot over (U)}]T from the host state filter 52.1. The position, velocity and acceleration components of the host vehicle 12 are calculated as follows:
θk=θk−1+ωk−1·T+{dot over (ω)}k−1·T2/2 (122)
xk=xk−1+cos θk·(Uk−1·T+{dot over (U)}k−1·T2/2) (123)
yk=yk−1+sin θk·(Uk−1·T+{dot over (U)}k−1·T2/2) (124)
{dot over (x)}k={dot over (U)}k·cos θk (125)
{umlaut over (x)}k={dot over (U)}k·cos θk+Uk·ωk·sin θk (126)
{dot over (y)}k=Uk·sin θk (127)
ÿk={dot over (U)}k·sin θk+Uk·ωk·cos θk (128)
Then from the clothoid model and equation (8):
After finding xk, yk and their derivatives using equations (122) through (128), equations (129) through (131) are used with curve fitting to solve for the road curvature parameters C0 and C1, wherein a window of sampled data points—e.g. in one embodiment, about 12 sample points—is used to improve the smoothness of the associated curve used for curve fitting.
A ninth embodiment of the curvature estimator 70 is a modification of the eighth embodiment, wherein the length L of the sliding window is adaptive. For example, the window length can be adapted to be responsive to the road curvature parameters C0 and C1, for example, in accordance with the following rule:
Lk=min{max{25−floor[(|C0m(k)|+33 |C1m(k)|)×7000],2},Lk−1+1} (132)
This rule provides for a larger window length L—e.g. as large as 25 samples—when both C0 and C1 are relatively small, for example, corresponding to a straight section of road. The window length L becomes smaller—e.g. as small as 2 samples—when either C0 or C1 is large, corresponding to a turn or transition of road. Furthermore, the window length L can be sharply reduced to account for a sudden vehicle maneuver, with limitation on the subsequent increase in window length L to one sample per step so that the previous samples so not adversely affect the output from the curvature estimator 70 when maneuver ends.
The above-described first through ninth embodiments of the road curvature estimation subsystem 42 are based upon the clothoid model of road curvature, wherein the road curvature C is characterized as varying linearly with respect to path length along the road, and wherein different types of roads (e.g. straight, circular or generally curved) are represented by different values of the clothoid road curvature parameters C0 and C1. The clothoid model reduces to a simpler form for straight (C0=C1=0) and circular (C1=0) road segments. Different embodiments of the road curvature estimation subsystem 42 can be used under different conditions. For example, in one embodiment of the predictive collision sensing system 10 in which the sampling rate of the yaw rate sensor 16 is relatively low, the seventh embodiment of the curvature estimator 70 is used when the longitudinal speed U of the host vehicle 12 is greater than a threshold, e.g. about 11 meters/second, and the ninth embodiment of the curvature estimator 70 is used at greater velocities. The ratio of the mean prediction error to the mean prediction distance can be used to compare and evaluate the various embodiments of the curvature estimator 70.
Generally, high-speed roads can be modeled as a collection of different types of interconnected road segments, each of which is represented by a different road model. For example,
Referring to
Accordingly, the roadway is represented by a multiple model system 82, wherein a particular road segment is characterized by one of a finite number r of models. For example, in
More particularly, the straight road model is characterized by:
C=0=>(C0=0 and C1=0) (133)
The circular arc or quadratic road model is characterized by:
C=C0=>(C1=0) (135)
where Tk is sampling period, Ûk is estimated host speed, and q=0.0005.
The clothoid road model is characterized by:
C=C0+C1·l (137)
where Tk is sampling period, Ûk is estimated host speed, and q=0.0025.
The multiple model system 82 has both continuous (noise) uncertainties as well as discrete (“model” or “mode”) uncertainties, and is thereby referred to as a hybrid system. The multiple model system 82 is assumed to be characterized by a base state model, a modal state, and a mode jump process.
The base state model is assumed to be characterized as follows:
x(k)=F[M(k)]x(k−1)+v[k−1, M(k)] (139)
z(k)=H[M(k)]x(k)+w[k, M(k)] (140)
where M(k) denotes the mode at time k in effect during the sampling period ending at k.
The modal state, or mode, is assumed to be one of r possible modes:
M(k)ε{Mj}rj=1 (141)
wherein the structure of the system and/or the statistics of the associated noise components can be different for different modes, as follows:
F[Mj]=Fj (142)
v(k−1, Mj)˜N(uj, Qj) (143)
The mode jump process governs the transition from one mode to another is assumed to be characterized by a Markov chain with known transition probabilities, as follows:
P{M(k)=Mj|M(k−1)=Mi}=pij (144)
The curvature estimators 70.1, 70.2, 70.3 operate in parallel, and the output therefrom is operatively coupled to a curvature processor 84—which, for example, can be embodied in the signal processor 26—which generates a single estimate of road curvature and associated covariance, in accordance with an interacting multiple model algorithm 2400 (IMM). Generally, the interacting multiple model algorithm 2400 is useful to track either or both maneuvering and non-maneuvering targets with moderate computational complexity, wherein a maneuver is modeled as a switching of the target state model governed by an underlying Markov chain. Different state models can have different structures, and the statistics of the associated process noises of different state models can be different. The interacting multiple model algorithm 2400 performs similar to the exact Bayesian filter, but requires substantially less computational power. Each model has a corresponding filter that is processed by the curvature processor 84 in accordance with the interacting multiple model algorithm 2400.
Referring to
Then, in step (2404), the mode-conditioned state is propagated for each model, i.e. j=1, . . . , r, according to a Kalman filter matched to the jth mode, Mj(k), so as to provide the state xj(k|k−1) and covariance Pj(k|k−1) at time k.
Then, in step (2406), the propagated mode-conditioned state estimates and covariances for each of the modes are combined to give:
Then, in step (2408), for each of the r parallel filters, the state estimates and covariances are calculated, conditioned on a mode being in effect and conditioned on the corresponding mode likelihood function. The Kalman filter matched to the jth mode, Mj(k), uses measurement z(k) to provide the state xj(k|k) and covariance Pj(k|k), whereby the likelihood function corresponding to filter j is given by:
Λj(k)=N[z(k); zj(k|k−1), Sj(k)] (149)
Then, in step (2410), the model probability μj(k) is updated for each mode, i.e. j=1, . . . , r, in parallel so as to provide the likelihood that each model is correct. The mixing probabilities are calculated for all combinations of (i,j) for (i,j=1, . . . r), as follows:
The mode probabilities are then updated for each mode, i.e. j=1, . . . , r, as follows:
Then, in step (2412), the overall state estimate and covariance—the output of the interacting multiple model algorithm 2400—is calculated by combining the mode-conditioned state estimates and covariances as follows:
The measurement and its noise covariance used as input to the interacting multiple model algorithm 2400 curvature filter, for each of the three associated mode-conditioned Kalman filters, are given by:
RCk=JCXPk|kJCXT (157)
where Pk|k is the state error covariance matrix from the host filter, and JCX is the Jacobian matrix given by:
The Markov model is implemented by assuming that at each scan time there is a probability pij that the target will make the transition from model state i to state j. These probabilities are assumed to be known a priori and can be expressed in a model probability transition matrix, e.g. as follows:
Referring to
[Xc, Yc, CR]=fM(x, y) (160)
wherein (x, y) are the world absolute coordinates of the vehicle position measure Zkg from the vehicle navigation system 86; Xc and Yc are vectors containing the coordinates of the center of the road closest to [x, y], and CR is an array containing the curvature parameters corresponding to the road center point coordinates in the vectors Xc and Yc. Accordingly, [Xc(i), Yc(i)] represents a point on the center of the road, and CR(i)=[C0, C1] represents the curvature parameters of the road at that point. The error covariance of Xc and Yc is Rmg, and the error covariance of CR is given by RCg, representing the noise or error in the associated map data. The map system 88 can either store both the road position coordinates and associated curvature parameter information, or could calculate the curvature parameters from the stored position coordinates of the road centers as the information is required. The world absolute coordinates can be transformed to host vehicle absolute coordinates by a combination of translation and rotation, given the position of the host vehicle 12 from the vehicle navigation system 86, and the heading of the host vehicle 12 based upon either information from the host filter 52.1, or from the vehicle navigation system 86, or both.
More particularly, given the vehicle position measure Zkg, in one embodiment the map system 88 uses an associated map database to determine the road that the host vehicle 12 is most likely on by finding roads in the associated map database that satisfy the following selection criteria:
(Zkg−Zr)·Rg·(Zkg−Zr)T≦T (161)
where T is a selection threshold and Zr is a point on the road that is closest to the vehicle position measure Zkg. If more than one road satisfies the selection criteria, then the most likely road is selected by comparing the curvature
estimated by another road curvature estimation system 42 with the curvature CR of each prospective road from the map database of the map system 88, and the road from the map database, for which the curvature CR of the point nearest to the vehicle position measure Zkg is closest to estimated curvature ĈkC from the other curvature estimation system 42, is selected as the most likely road from the map database, and the associated curvature at the closest point is then given as Ck|kg=CR.
Referring to
and an associated covariance Pk|kC
Referring to
of the target vehicle 36—as measured by the radar sensor 14—can be used to determine an estimate of road curvature of the roadway 34 based upon the premise that under normal driving conditions, the target vehicle 36 is assumed to follow the roadway 34. The dynamics of the tth target are assumed to be given by the following constant-acceleration kinematic equations, which are embodied in an extended Kalman filter 90 of an auxiliary filter 90.1.
xk+1at=Fkat·xkat+wkat, wkat˜N(0,Qkat) (162)
zkat=hkat(xkat)+vkat, vkat˜N(0,Rkat) (163)
where
wherein T is the sampling period, the superscript (•)t is used to designate a particular target, and the superscript (•)a is used to indicate that the filter is auxiliary.
The associated measurement function hkat is non-linear, and give n by:
The extended Kalman filter 90 provides for linearization of the non-linear measurement function hkat using an associated Jacobian matrix, and provides for an estimate of the target state {circumflex over (x)}k|kat and its associated error covariance Pk|kat, which are transformed by an associated measurement processor 92 so as to provide the measurement input zkC
RkC
The system equations of the associated Kalman filter 94 of the associated curvature filter 94.1 are given by:
xk+1C
zkC
where
and the associated measurement matrix HkC
The auxiliary filter 92.1 and the curvature filter 94.1 of the road curvature estimation subsystem 42″ operate in accordance with the Appendix.
Referring to
From the first road curvature estimation subsystem 42.1, the curvature estimate and associated error covariance at a distance l along the roadway 34 from the current location are given respectively by:
C(l)=C0+C1·l (172)
RC(l)=RC0+RC1·l (173)
where RC0, RC1 and RC are the error covariances of C0, C1 and C respectively.
From the second road curvature estimation subsystem 42.2, the corresponding curvature estimate is given by Cg(l) and the corresponding error covariance is given by Rg, wherein Rg is generally a constant scalar.
The curvature estimates and associated error covariances are combined by the road curvature fusion subsystem 96, for example, as follows:
Cf(l)=G·C(l)+Gg·Cg(l) (174)
Rf(l)=(RC(l)−1+Rg−1)−1 (175)
where G and Gg are weights given by:
The fused curvature Cf(l) and associated error covariance Rf(l) can be used by other processes, for example for improving the estimation of the locations of the host 12 or target 36 vehicles, or for collision prediction. For example, in the embodiment of
While specific embodiments have been described in detail in the foregoing detailed description and illustrated in the accompanying drawings, those with ordinary skill in the art will appreciate that various modifications and alternatives to those details could be developed in light of the overall teachings of the disclosure. Accordingly, the particular arrangements disclosed are meant to be illustrative only and not limiting as to the scope of the invention, which is to be given the full breadth of the any claims derivable from the description herein, and any and all equivalents thereof.
A Kalman filter is used to estimate, from a set of noisy measurements, the state and associated covariance of a dynamic system subject to noise.
The system dynamics are defined by:
xk+1+Fk·xk+wk, wk˜N(0, Qk) (A-1)
where xk is the system state vector, Fk is the system matrix and wk an associated vector of noise variables corresponding to each state variable, each noise variable having a mean value of zero, and a variance given by the corresponding element of the associated variance vector, Qk.
The dynamics of the associated system measurements are given by:
zk=Hk·xk+vk, vk˜N(0, Rk) (A-2)
where zk is the system measurement vector, Hk is the measurement matrix and vk an associated vector of noise variables corresponding to each measurement variable, each noise variable having a mean value of zero, and a variance given by the corresponding element of the associated variance vector, Rk. The values of the elements of the associated covariance matrix Rk can be determined a priori from analysis of the representative measurements of the associated system for associated representative sets of operating conditions. The values of the elements of the associated covariance matrix Qk account for modeling errors. Generally, the associated matrices Fk, Qk, Hk, Rk can vary over time.
Given a measurement zk at time k, and initial values of the state xxk−1|k−1 and associated covariance Pk−1|k−1 at time k−1, the Kalman filter is used to estimate the associated state xk|k and associated covariance Pk|k at time k.
The first step in the filtering process is to calculate estimates of the state xk|k−1 and associated covariance Pk|k−1 at time k based upon estimates at time k−1, as follows:
xk|k−1=Fk·xk−1|k−1 (A-3)
PR|k−1=Fk·Pk−1|k−1·FkT+Qk (A-4)
The next step is to predict the measurement {circumflex over (z)}k and associated covariance matrix Sk at time k, as follows:
{circumflex over (z)}k=Hk·xk|k−1 (A-5)
Sk=cov({circumflex over (z)}k)=Hk·Pk|k−1·HkT+Rk (A-6)
The next step is to calculate a gain matrix Gk used for updating the state vector xk|k and associated covariance matrix Pk|k, as follows:
Gk=Pk|k−1·HkT·Sk−1 (A-7)
Finally, the state vector xk|k and associated covariance matrix Pk|k are estimated at time k, responsive to the associated measurement zk, as follows:
xk|k=xk|k−1+Gk·(zk−{circumflex over (z)}k) (A-8)
Pk|k=Pk|k−1−Gk·Sk·GkT (A-9)
The instant application is a continuation-in-part of U.S. application Ser. No. 10/620,749 filed on 15 Jul. 2003, now U.S. Pat. No. 7,034,742 which claims the benefit of prior U.S. Provisional Application Ser. No. 60/396,211 filed on Jul. 15, 2002. The instant application also claims the benefit of prior U.S. Provisional Application Ser. No. 60/532,344 filed on Dec. 24, 2003. The above-identified applications are incorporated herein by reference in their entirety.
Number | Name | Date | Kind |
---|---|---|---|
2709804 | Chance et al. | May 1955 | A |
3177485 | Taylor, Jr. | Apr 1965 | A |
3603994 | Williams et al. | Sep 1971 | A |
3699573 | Andrews et al. | Oct 1972 | A |
3725918 | Fleischer et al. | Apr 1973 | A |
3869601 | Metoalf | Mar 1975 | A |
3971018 | Isbister et al. | Jul 1976 | A |
4623966 | O'Sullivan | Nov 1986 | A |
5051751 | Gray | Sep 1991 | A |
5138321 | Hammer | Aug 1992 | A |
5170440 | Cox | Dec 1992 | A |
5202691 | Hicks | Apr 1993 | A |
5307289 | Harris | Apr 1994 | A |
5314037 | Shaw et al. | May 1994 | A |
5343206 | Ansaldi et al. | Aug 1994 | A |
5402129 | Gellner et al. | Mar 1995 | A |
5406289 | Barker et al. | Apr 1995 | A |
5471214 | Faibish et al. | Nov 1995 | A |
5530651 | Uemura et al. | Jun 1996 | A |
5537119 | Poore, Jr. | Jul 1996 | A |
5539397 | Asanuma et al. | Jul 1996 | A |
5587929 | League et al. | Dec 1996 | A |
5594414 | Namngani | Jan 1997 | A |
5598164 | Reppas et al. | Jan 1997 | A |
5627768 | Uhlmann et al. | May 1997 | A |
5631639 | Hibino et al. | May 1997 | A |
5633642 | Hoss et al. | May 1997 | A |
5657251 | Fiala | Aug 1997 | A |
5668739 | League et al. | Sep 1997 | A |
5684473 | Hibino et al. | Nov 1997 | A |
5689264 | Ishikawa et al. | Nov 1997 | A |
5703593 | Campbell et al. | Dec 1997 | A |
5710565 | Shirai et al. | Jan 1998 | A |
5751211 | Shirai et al. | May 1998 | A |
5926126 | Engelman | Jul 1999 | A |
5948043 | Mathis | Sep 1999 | A |
5955967 | Yamada | Sep 1999 | A |
5959552 | Cho | Sep 1999 | A |
5959574 | Poore, Jr. | Sep 1999 | A |
5983161 | Lemelson et al. | Nov 1999 | A |
6070121 | Matsuda | May 2000 | A |
6085151 | Farmer et al. | Jul 2000 | A |
6134509 | Furusho et al. | Oct 2000 | A |
6141617 | Matsuda et al. | Oct 2000 | A |
6161071 | Shuman et al. | Dec 2000 | A |
6198426 | Tamatsu et al. | Mar 2001 | B1 |
6226389 | Lemelson et al. | May 2001 | B1 |
6265991 | Nishiwaki et al. | Jul 2001 | B1 |
6268825 | Okada | Jul 2001 | B1 |
6275231 | Obradovich | Aug 2001 | B1 |
6275773 | Lemelson et al. | Aug 2001 | B1 |
6282478 | Akita | Aug 2001 | B1 |
6292752 | Franke et al. | Sep 2001 | B1 |
6300865 | Fechner et al. | Oct 2001 | B1 |
6343253 | Matsuura et al. | Jan 2002 | B1 |
6370475 | Breed et al. | Apr 2002 | B1 |
6374184 | Zahm et al. | Apr 2002 | B1 |
6393376 | Andreas | May 2002 | B1 |
6405132 | Breed et al. | Jun 2002 | B1 |
6420997 | Cong | Jul 2002 | B1 |
6424904 | Takahashi et al. | Jul 2002 | B1 |
6470272 | Cong et al. | Oct 2002 | B2 |
6542111 | Wilson | Apr 2003 | B1 |
6567039 | Shirai et al. | May 2003 | B2 |
6643588 | Ibrahim | Nov 2003 | B1 |
6670910 | Delcheccolo et al. | Dec 2003 | B2 |
6675090 | Matsuura | Jan 2004 | B2 |
6675094 | Russell et al. | Jan 2004 | B2 |
6718259 | Khosla | Apr 2004 | B1 |
6751547 | Khosla | Jun 2004 | B2 |
6763318 | Winter et al. | Jul 2004 | B1 |
6823241 | Shirato et al. | Nov 2004 | B2 |
20020022927 | Lemelson et al. | Feb 2002 | A1 |
20020042668 | Shirato et al. | Apr 2002 | A1 |
20020044080 | Shirai et al. | Apr 2002 | A1 |
20020049539 | Russell et al. | Apr 2002 | A1 |
20030100992 | Khosla | May 2003 | A1 |
20030195703 | Ibrahim | Oct 2003 | A1 |
20030218563 | Miyahara | Nov 2003 | A1 |
20040143416 | Hattori et al. | Jul 2004 | A1 |
20040183663 | Shimakage | Sep 2004 | A1 |
20050179580 | Cong et al. | Aug 2005 | A1 |
Number | Date | Country | |
---|---|---|---|
20050225477 A1 | Oct 2005 | US |
Number | Date | Country | |
---|---|---|---|
60532344 | Dec 2003 | US | |
60396211 | Jul 2002 | US |
Number | Date | Country | |
---|---|---|---|
Parent | 10620749 | Jul 2003 | US |
Child | 11022265 | US |