System and method for estimating a shooter range

Systems and methods to disambiguate positions of a shooter

Field of the Invention

The present invention relates to safety and technologies for law enforcement, and more specifically to systems for determining the origin and direction of movement of supersonic projectiles. The methods and systems are capable of examining and disambiguating positions of a shooter, even for long distances between the shooter and the sensor, and when no signal is received or only a weak signal is received from the sound of the muzzle.

Background of the invention

Systems and methods are known that can determine in general the direction and trajectory of supersonic projectiles, such as bullets and artillery shells, by measuring parameters associated with the shock wave generated by a projectile. One such system, described in US Pat. Number 5 241 518, includes at least three spatially separated sensors, each sensor incorporating three acoustic transducers arranged in a plane. The sensors generate signals in response to the shock wave, which are related to the azimuth and elevation angle with respect to the origin of the shock wave. With measurements of shock waves only, it is impossible to determine the distance between the sensor or sensors and the origin of the shock wave. The distance information is typically obtained from the flash or the burst.

The azimuth and the angle of elevation of a shooter in relation to the location of the sensor, are typically determined by measuring information of the Moment of Arrival (TOA) of the cannon mouth signal and the wave signal of shock, on each sensor. Each of the sensors finds the signals at different times, and generates a signal in response to the pressure of the muzzle and the shock wave. The signals from various sensors are processed, and a direction (azimuth and elevation) can be determined from the sensor or sensors to the origin of the barrel mouth and the shock wave, and therefore the projectile's trajectory.

Conventional systems use microphones, which can be relatively close (for example 1 meter apart) or widely separated (for example mounted on a vehicle, or carried by soldiers on a battlefield), and measure the pressure of the mouth of the cannon and shock wave omnidirectionally in their respective locations. However, unless the sensors have a relatively large separation and / or the path is within the antenna, the precision of synchronization necessary to obtain precise solutions only with shock wave is very high, and special techniques are required.

A large antenna can be a major disadvantage, for example in vehicle-mounted systems. In addition, systems with only marginal temporal resolution can produce ambiguous solutions, since the information of the moment of arrival of the shock wave in a given set of sensors is almost identical for two shooter positions in specular symmetry.

Conventional algorithms require at least 4 shockwave and cannon mouth detections, so that a 4x4 matrix can be inverted to map a flat wave over the TOA of the shockwave. Small errors in the determination of the TOA of the shock wave and of the muzzle can produce substantial errors in the scope of the estimates. In addition, conventional algorithms assume a constant bullet velocity along the trajectory of the bullet, which provides inaccurate range estimates for long-range shots produced from a distance of more than about 300 m.

US 5 930 202 discloses a method for determining an unambiguous trajectory of a projectile using a matrix of separate acoustic sensors, in accordance with the precautionary part of claim 1.

Summary of the Invention

The invention addresses the shortcomings of the prior art, providing a method characterized in accordance with claim 1.

Preferably, the method further includes executing a genetic algorithm for a predefined number of generations, with an initial population that includes a predetermined number of individuals, each individual being represented by a 4-tuple that includes the azimuth of the shooter, the elevation of the shooter, the failed azimuth and the failed elevation, and calculate residues for the individuals in each generation, including in the residuals a least squares adjustment of a combination of TDOA and rebufo shock wave signals. If the relationship between the solution with the minimum residue and its ambiguous alternative solution is greater than a predefined value, for example a value of at least 2, then the solution with the minimum calculated residue is designated as the trajectory of the disambiguous projectile. To prevent spurious signals from being interpreted, such as shockwave waveforms, a projectile path can be eliminated as false if the acoustic energy of the measured shockwave waveform has less than a threshold value above a predetermined frequency band, for example frequencies between approximately 700 Hz and 10 kHz. Alternatively or additionally, the

5 trajectory of a projectile and be considered false, if an estimated interval in which a waveform of the shock wave has a positive value, is less than a minimum time or greater than a maximum time, for example less than about 70 μs or greater than approximately 300 μs.

In advantageous embodiments, the total energy can be determined by integrating the energy emitted on the window, preferably neglecting parts of the detected signal, caused by echoes of the shock wave. 10 Advantageously, the maximum signal value can be determined in the window that produces the maximum total energy and, if the maximum signal value is greater than the total energy measured in the window by a preferred proportional factor, the maximum signal value can be identified as being related to the signal from the muzzle.

Embodiments of the invention may include one or more of the following characteristics. The distribution of synchronization error, of the antenna and / or the acoustic sensors, can be related to gain variations,

fifteen sampling variations and variations of sensor location, of antenna sensors. The level of confidence for disambiguation depends on the size of the antenna, so that smaller antennas require greater measurement accuracy. If there are two ambiguous solutions, the disambiguous trajectory of the projectile is selected based on a ratio of the waste to two ambiguous solutions.

In other advantageous embodiments, differences between arrival times (TDOA) for sensor pairs may

twenty determined by designating as the reference sensor a sensor that first finds the shock wave, and activating a first safety of a synchronization circuit when the amplitude, for example the initial part of the shock wave only signal, in the sensor reference, crosses a threshold value. The first insurance activates start counters for each of the other sensors, the counter running on each of the other sensors until the corresponding sensor finds the shock wave. When one of the other sensors finds, for example, the

25 Initial part of the shock wave only signal activates a second safe for such a sensor, which stops the initial counter for such a sensor. The TDOA values ​​are then recorded for the other sensors, in relation to the reference sensor.

Other features and advantages of the present invention will be apparent from the following description of the preferred embodiments, and from the claims.

30 Brief description of the drawings

These and other features and advantages of the invention will be more fully understood by means of the following illustrative description, referring to the attached drawings, in which elements are labeled with reference numbers, and which may not be to scale.

Figure 1 schematically shows a cross-sectional view of a Mach cone crossing an antenna;

35 Figure 2 schematically shows a sensor array with 7 omnidirectional acoustic sensors;

Figure 3 schematically shows the ambiguity inherent in determining the path only with a shock wave;

Figure 4 schematically shows a probability density for the time difference in the arrival measurements, for the determination of the curvature of the Mach cone;

40 Figure 5 schematically shows the probability of correct disambiguation, between shooter trajectories;

Figure 6 shows a schematic diagram of a correlation process;

Figure 7 is a process flow of a genetic algorithm, used to correctly disambiguate between shooter trajectories;

Figure 8 is a process flow to discriminate against non-shock wave signals;

Four. Five Figure 9 is a schematic diagram of a shock wave arrival time (TOA) model;

Figure 10 shows a schematic diagram of process flow, for scope estimation; and Figure 11 shows a schematic diagram of process flow of a genetic algorithm for scope estimation.

Detailed description of certain illustrated embodiments

As described above in summary form, the invention provides in various embodiments, methods and systems for estimating the range of a shooter and disambiguation of projectile trajectories. These systems and methods are especially useful and advantageous when an insufficient number of parameters required for a precise solution is detected, or when such parameters cannot be detected reliably.

The trajectories of supersonic projectiles are estimated exclusively from moments of arrival of the shock wave of the projectile, measured by several sensors close to each other, distributed by a "small" measuring volume, called an antenna. A measurement volume is considered small if the sensor separation is 2 meters or less. Once the trajectory of the projectile is identified, the position of the shooter is known, except for the return distance along the trajectory. This distance can be found if the antenna also gets the moment of arrival of the rebufo sound. However, the rebufo is not always detectable, so an accurate solution from shock wave only is essential to determine the trajectory.

Next with reference to Figure 1, the surface of the shock wave is considered to be a conical expansion surface, which has its axis coinciding with the trajectory of the bullet. The surface of the shock wave is also called the Mach cone. To obtain the solution only by shockwave, three properties have to be determined from the arrival times measured in five or more antenna sensors, namely the angle of arrival, the radius of curvature and the spatial gradient of the radius of curvature of the conical expansion surface.

The angle of arrival of the conical surface generator that first reaches the antenna determines two possible relative angles (often called 'ambiguous' angles) of the bullet's trajectory, relative to the angle of arrival at the antenna. The 'ambiguous' angles will be described in greater detail below, referring to Figure 3. The radius of curvature of the conical surface in the antenna determines both the distance and the direction of the path. The gradient of the radius of curvature along the surface generator path determines in which direction the bullet moves, thereby eliminating the 'ambiguity' between two possible directions. Determining precisely these three properties of the shock wave and deciding correctly between two possible 'ambiguous' trajectories requires very precise measurements. For example, random errors should not be greater than about 1 μs to correctly decide between the two alternate angles of orientation of the shooter.

The necessary accuracy can be estimated by considering the propagation characteristic of the shock wave described in Figure 1. Next also referring to Figure 2, an antenna 20 includes N sensors (N = 7) capable of determining the arrival times of a conical shock wave that advances. Since it is expected that the incoming trajectories of the bullet are essentially originated anywhere, the antenna elements 23 to 28 may advantageously be distributed uniformly at locations C (Cxj, Cyj, Czj) on a spherical surface, with an element 22 located in the center at (Cx0, Cy0, Cz0), so that there is a uniform sensor opening regardless of the angle of arrival. The time in which the first sensor, designated as the reference sensor, detects the conical surface in advance, is denoted as to. The other sensors detect the conical surface in advance in subsequent defeated moments like you. The sound propagation distances in the direction of the conical feed surface are obtained by multiplying each of the time differences by the local velocity of the sound c, that is di = c • (ti - to). If there are no measurement errors, then the conical surface that passes through the reference sensor is also determined by the other (N - 1) sensors, with the three-dimensional coordinates of the N points, ideally determining all the parameters of the cone of the shock wave. However, as mentioned above, errors in arrival time measurements and sensor coordinates can result in erroneous parameters for the shock wave cone, and therefore also for the projectile's trajectory. In the following, the necessary details will be described, in the difference between arrival times, to make correct decisions around the two trajectory angles, otherwise ambiguous.

Advantageously, the system incorporates features that ensure that it will not confuse non-ballistic signals such as vehicle noise, vibration, wind noise and EMI (electromagnetic interference), with a shooter. For example, the sensor mast may be mounted on a vehicle (not shown) with elastomer sleeves in coupling joints, to prevent rattling. The sensors can be attached to the ends of the spins with elastomeric couplings, which have low frequency resonances at approximately 1 Hz, to isolate them from spin vibration. The sensor spins can be attached to a common core that contains analog electronics, which can also be attached to the sensor mast with elastic brackets to isolate it from mast vibrations.

In addition, the following decision algorithm can be used to filter signals that lack the distinctive characteristics typically found in signals derived from shock wave. All values ​​are parameterized, that is, they are relative and can be adjusted externally. The listed values ​​are provided by way of illustration only.

Referring to Figure 8, a process 800 determines whether a detected signal originates from a shock wave. Process 800 begins at step 802, and at step 804 it checks if the signal is an event.

5 strong enough to be considered a shock wave, for example if the maximum signal value exceeds a specific parameterized threshold, for example 500. If this is the case, process 800 continues with step 806 and checks if there is an abrupt transition from zero to the maximum signal value, ensuring that the transition to this maximum value is not preceded by another signal having a significant magnitude, for example 1/16 of the maximum signal value.

10 If this is the case, process 800 continues with step 808 and verifies if the time between the minimum and the maximum of the shock wave has a sufficiently large value, for example 200-400 μs. If this is the case, process 800 goes to step 810 and checks if the quantities are in close proximity, for example within 35%. If this is the case, process 800 continues in step 812 and verifies whether the maximum transient pressure between the minimum and zero peak signal is abrupt, using essentially the same criteria as in step 806. If this is the case, he

fifteen process 800 continues with step 814 and verifies whether the times between the maximum signal value and the zero crossing are comparable, and between the zero crossing and the minimum signal value, for example within approximately 180 μs. If all stages produce an affirmative response, process 800 decides that the signal may be a shock wave and the signal is processed, step 816. Conversely, if the response to one of the 6 decision stages is negative, the detected signal has not originated from a shock wave, step 818.

twenty Referring again to Figure 1, it is assumed that the trajectory of the projectile coincides with the x-axis. The angle of Mach is given by 8 = arcsen (1 / M), where N is the number of Mach, defined as the velocity of projectile V divided by the velocity of sound c. L refers to the characteristic length of the antenna. The radii of curvature of the cone at the two ends of the antenna 20 are r1 and r2. The side view of the left half of the drawing shows how the curvature r1 is measured. The distance d is equal to ad = r1 • cos (<). The angle <is defined by sin (<) =

25 L / 2 r1, so that for small angles of <F ~ L / 2 r1 is obtained. The measure of the temporal difference of the curvature between the points of the antenna surface that bisect the conical surface with radius r1, is equal to dt1 = / d / c = (r1 - d) / c ~ r1 <2 / 2c = L2 / (8 • r1 • c). The measure of temporal difference of the curvature in r2 = r1 - L • sin (8), is given by the same expression substituting r2 for r1. Therefore, dt2 = dt1 + L3 sin (8) / 8 r12 c.

Assuming non-biased measurement errors, that is, assuming that the temporal differences of measurement dt1 and dt2 30 are randomly distributed values ​​that have different averages dt1 and dt2 but the same standard deviation

or statistically determined, the averaged measurement values ​​at the two ends of the matrix correctly determine the local curvature at that point. Example distributions of measured values ​​for the temporal differences dt1 and dt2 are shown in Figure 4.

The sample measurement taken at end 2 is shown as X. The radius of curvature at end 2 (radius r2)

35 It is less than one end 1 (radius r1). Therefore, all measurements taken at end 1 that have values ​​greater than X will result in the correct decision that the curvature at end 1 is greater than at end 2. The probability of making the right decision when the measure at end 2 it is equal to X is given by:

40 Integrating into x and performing variable substitution results in the following probability of making the right decision:

Referring to Figure 5 below, the probability of a correct decision, or the level of confidence for disambiguation, is represented for two example antenna sizes, L = 1m and L = 2m, versus the nearest point Approach (CPA) r between the trajectory of the projectile and the antenna

twenty. It is assumed that the speed of sound is c = 340 m / s. It is evident that a larger antenna has a significantly expanded range, for non-ambiguous solutions only with a shock wave. For large CPA values, the difference in curvature at the two ends of the antenna (r1 and r2) is too small to be distinguishable, so that the probability of a correct decision is close to 50%, that is to say ambiguity It is total. Therefore, the level of confidence depends on the size of the antenna, that is to say its diameter or its spatial extent.

As mentioned above, errors are generated from synchronization errors and uncertainty in the sensor coordinates. Uncertainty in the coordinates of the sensor contributes to systematic errors, which are an extremely variable function of the angle of arrival of the shock wave. However, for random arrival angles the sensor coordinate errors appear as random errors for the time difference.

Synchronization errors also arise from variations in both gain and signal strength, between one channel and another. Arrival times are obtained when the sensor outputs are increased to a predetermined threshold value V0. The synchronization error dt caused by a variation of gain dg, depends

of the temporary rate of increase of tension for the channel, having dg V0.

dt = g dV / dt

Synchronization errors also occur when the signal strength varies over the aperture. For an opening of length L and a cylindrical sound source at distance r, the maximum variation of the signal level through the opening is equal to p0 (L / 2r), where p0 is the sound pressure at the center of opening. The previous synchronization error equation also applies to this type of error, with the expression L / 2r replacing the relative gain variation dg / g. Amplitude errors are not random between the sensors, but vary uniformly from a maximum in the full opening to zero in the center. For ranges greater than 10 m, for an opening of 1 m the maximum amplitude factor is less than 0.05, which is less than the channel gain variation parameter of 0.2, so that due effects can be ignored to amplitude errors. Conversely, as described above, for scopes smaller than about 10 m the radius of the Mach cone is small enough with respect to the opening length of 1 m so that the measurement errors are not significant.

With realistic estimates for synchronization errors caused by sensor uncertainty, assuming that the magnitudes of the error vectors are statistically independent and uniformly distributed between 0 and 1 mm, and that the error angles are statistically independent, the standard deviation of uniformly distributed random time differences will be equal to

The standard deviation of random time sampling errors, distributed in a binomial manner, for a 1 MHz system sampling, is equal to 0.25 μs. It is estimated that synchronization errors due to gain variations are approximately 0.75 μs, for an exemplary system with a channel bandwidth of approximately 18 kHz, corresponding to a voltage variation of approximately 0.02 V / μs. The acoustic sensors used for each matrix are chosen with sensitivities within ± 1.5 dB. Therefore, the relative gain variations of the channel are distributed approximately evenly between 0.84 and 1.19, so that the standard deviation of the relative gain is

approximately equal to

The threshold voltage is V0 = 0.15 V, which results in a standard deviation of synchronization errors of approximately 0.75 μs.

The total measurement synchronization errors are estimated assuming that channel gain variations, sampling variations and sensor location variations are all statistically standard synchronization error can be estimated as

It is complicated and expensive to achieve such precision with analog to digital conversion, because high sampling frequencies followed by interpolation are required. In the system disclosed, two different circuits are used to accurately measure the difference between arrival times (TDOA).

In one embodiment, the exemplary system uses an analog circuit of difference between arrival times (TDOA), using 1 MHz clocks on each channel. The clocks are triggered when the sensor signal exceeds a threshold signal level in the reference sensor, which was defined above as the sensor that first encounters the shock wave. As discussed above, a clock frequency of 1 MHz is sufficient to eliminate the importance of temporary sampling errors in practice. The system operates in analog mode, based on the detection of threshold levels, with digital logic performing the following functions:

5 1. A first safe is activated when the signal amplitude of the channel, in the reference sensor that first encounters the shock wave, crosses a threshold value.

<dl><dt>2. </dt><dd>The first insurance activates start counters for each channel, which increase your account by one in each clock cycle. The processor is alerted.</dd></dl>

<dl><dt>3. </dt><dd>The counter on each channel works until the corresponding sensor finds the shock wave. This</dd></dl>

10 activates a second safe on the channel, which stops the account on that channel. If a second safe is not activated, the corresponding counter operates up to a higher limit value.

<dl><dt>4. </dt><dd>The final number of computations in each counter is registered in a TDOA digital register. </dd></dl>

<dl><dt>5. </dt><dd>The processor reads the TDOA record. </dd></dl>

<dl><dt>6. </dt><dd>The processor resets the counters to receive the next shock wave. </dd></dl>

fifteen In another embodiment, the correlation is calculated for each channel with all other channels, for a time segment centered on the instant of detection of the TDOA by physical equipment. The correlation of two functions, denoted Corr (g, h), is defined by

The correlation is a function of t, which is called a "offset." Therefore it is located in the time domain and 20 has the following property:

where gyh are real functions of time. G (f) is the Fourier transform of g (t), and H (f) is the Fourier transform of h (t).

The total power in a signal is:

The arrival time signal has a finite length, so it is only necessary to perform the integration (or sum, for discrete data) over a finite time interval, centered around the arrival time; The length of the data in one or both channels can be extended by zero padding, so that the duration of the two signals matches, as is known in the art.

30 In the following discussion, integral functions of continuous functions are used for simplicity, although the actual data are digitized and discrete values. Those skilled in the art will be able to easily replace the integrals with a sum.

Next, referring to Figure 6, in a process 60 the temporal data of the shock wave signal gi (t), gj (t) are acquired on each channel i, j, step 601, 602, and recorded according to weather. In the stages

35 603, 604, the total signal power on a channel i is calculated for the subsequent normalization of the correlation, such as

The Fourier transform Gi (f) of the temporal data of the shock wave signal gi (t) is calculated for channel i, and the conjugate Gi (-f), step 605 is obtained. Similarly, the transform of Fourier Gi (f) of the temporal data of the shock wave signal gi (t) is calculated for all other channels j, step 606. Next, the cross correlation Gi (-f) is obtained · Gj (f) for each pair of channels (i, j), step 608, which is a fi, j (t) function of the "delay" t. The TDOA for each pair of channels, is the time tmax at which f (t) has its maximum value, step 610. The correlation between channels i and j can be defined as

The residue for channel i is calculated by calculating the average value of a sensor i over all sensors j:

as indicated in step 612. The TDOA and correlations for that channel with the best (ie the smallest) overall residue are then selected as the "best" solution, step 614.

As mentioned above, channel data is typically sampled at discrete time intervals, with a predefined sampling rate, for example 41,666.66 samples / second. This

fifteen corresponds to a width of 24 μs, which reflects the temporal resolution for the received signal. The correlation process is carried out with a temporal resolution that is improved by a factor of 8, up to 3 μs, by taking 333 333 samples / second.

Once the various differences between arrival times (TDOA) between the sensors have been determined, based on shock wave signals only, the azimuth and elevation of the shooter can be determined, as well as the trajectory of

twenty the bullet. The position of the shooter, that is the distance of the shooter from the sensor array, can be determined if the signal from the rebufo is also known.

In a Cartesian coordinate system centered in the center of the matrix, ie {(Cx0, Cy0, Cz0) = (0, 0, 0)}, the TOA arrival time of the shock wave at a given sensor (Cxj , Cyj, Cyj) (see Figure 2), is given by:

Vx

represents the supersonic speed of the bullet

where c is the speed of

V

Vy

=

Vz

sound and N the number of Mach. 1 represents the 'failed angle' between the position of the shooter and the trajectory of the bullet, which includes both the azimuth and elevation angles. A direct impact would correspond to 1 = 0. The angle of

sen ()

1

Mach 8 is defined by

=

.

M

As mentioned above and indicated in Figure 3, for a given shooter position and a bullet path 30, there are other shooter positions and bullet trajectory, for which the TOA of the shock wave in a set of sensors is almost identical. In fact, the two ambiguous solutions are identical if, in a simplified model, it is assumed that the shock wave propagates through the sensor array as a flat wave. If the resolution of the TDOA is high enough to resolve the curvature of the shock wave, then the two almost identical solutions can be disambiguated. The essential ambiguity of TDOA solutions only with a shock wave is indicated in Figure 3.

Assuming precise enough TDOA measurements, the true solution for the shooter's position and the trajectory of the bullet can be obtained by calculating the shooter / trajectory combination that minimizes the residual mean square value (RMS) of the TDOAs measured and calculated shock wave:

10 where the sum is made on all the sensors.

One approach to solve this problem is the Levenberg-Marquardt L1 algorithm, described in detail in US Pat. 5 930 202. Most classic point-to-point algorithms use a deterministic procedure to approximate the optimal solution, starting at an assumed random solution, and specifying a search direction based on a previously specified transition rule, such as methods

fifteen direct that use an objective function and limit values, and gradient-based methods that use first and second order derivatives. However, these methods have disadvantages, for example that an optimal solution depends on the initial solution selected, and that the algorithm can be locked in a solution below the optimum level, such as a local minimum or where the surface of the function has a valley flat, so that additional iterations do not improve the result.

twenty It has been found that a global minimum of the direction of the shooter and the trajectory of the projectile can be calculated more quickly and disambiguated with greater reliability by using an evolutionary genetic algorithm (GA). GAs mimic natural evolutionary principles and apply them to search and optimization procedures.

A schematic flow diagram of a GA is shown in Figure 7. Instead of starting with a single solution hypothesis, a GA 70 process begins its search by initializing a random population of solutions,

25 step 71, and set a generation counter to zero, indicating the initial set of solutions, stage 72. Once a random population of solutions has been created, each is evaluated in the context of the non-linear programming problem, stage 73 , and an aptitude (relative merit) is assigned to each solution, step 74. The aptitude can be represented by the Euclidean distance / Tmin between a calculated solution and the measured solution.

30 Intuitively, an algorithm that has a small value of / Tmin is better. For example, when the GA is applied to disambiguate the solution for the direction of the shooter and the trajectory of the projectile, the GA as an example uses as an chromosome an initial population of 200 4-s, with each 4- containing the following values: [Azimuth Puller, Elevation Puller, Azimuth Failed, Elevation Failed]. 35 [Azimuth Puller, Elevation Puller] are defined by the angle (8 + 1), while [Azimuth Failed, Elevation Failed] are

defined by angle 1 (see figure 3). Since the rebufo is not used with the shock wave only approach described above, a nominal range of 100 meters is assumed between the sensor array and the shooter. The initial population is created by random 4-s selection, covering a reasonable and significant range of values

(all values ​​are given in degrees): 40 Shooting Azimuth = {0, ..., 360}, Shooting Elevation = {-10, ..., 30},

Azimuth Failed = {-20, ..., 20}, and

Failed Elevation = {-20, ..., 20}.

In step 75, it is checked whether a maximum number of iterations has been reached for the GA, which can be set, for example, 25. If the maximum number of iterations has been reached, process 70 stops at step 80, and 5 result can be accepted or evaluated additionally. In another case, step 76 verifies whether the predetermined proficiency criterion is met.

For example, the fitness criterion can be a calculated azimuth of &lt;15 ° and / or a ratio between the residues of two ambiguous solutions. If the fitness criterion is satisfied, process 70 stops at step 80; otherwise, a new population is created by crossing, stage 77, and mutation, stage 78, and the generation counter is

10 Increase by one, stage 79.

In each generation the "best" individual is allowed to survive without mutation, while the best 100 individuals estimated according to their ability also survive, but the following 100 individuals are used to create with the crossover / mutation operators listed in table one of pairs of these survivors.

The following crossover and mutation operators are used as an example to show process 70:

fifteen Table 1

<dl><dt>Name of the operator </dt><dd>Type of Operator Probability Description </dd></dl>

<dl><dt>Crossing - Azimuth </dt><dd>Crossing 0.5 Swap handle / path azimuth between two chromosomes </dd></dl>

<dl><dt>Crossing - Failed </dt><dd>Crossing 0.5 Exchange failed azimuth / elevation between two chromosomes </dd></dl>

<dl><dt>Mutation - Field </dt><dd>Mutation 0.3 Replace a given field (with a probability of 0.25 per field) with a new randomly selected value within range </dd></dl>

<dl><dt>Mutation - Incremental </dt><dd>Mutation 0.4 Induce small mutations in all fields of a chromosome (within: 2 ° for shooter information; within: 5 ° for failed information) </dd></dl>

<dl><dt>Mutation - By Alternation </dt><dd>Mutation 0.1 Change the solution in the ambiguous alternative solution </dd></dl>

<dl><dt>No - Mutation </dt><dd>Mutation 0.2 The chromosome remains intact </dd></dl>

The disambiguation is achieved and / or improved by carrying out a gradient search of the best solution and the corresponding alternative solution. For ambiguous solutions, waste and waste rates are calculated. If the calculated failed azimuth is &lt;15 °, which represents "next" shots, and if the ratio of the

twenty waste is&gt; 2, then the solution with the smallest residue is selected. In another case no selection is made, and the solution with the lowest residue is labeled as the "main" solution, the other solution being labeled as an "alternative" solution.

With detection only by shock wave, the GA algorithm produced a solution on a 1 GHz computer operating with the Linux operating system in 0.15 seconds, over a wide range of simulated shots. 97% 25 of the simulated shots are within 15 ° of the failed azimuth, and 86% of the simulated shots are within 5% of the failed azimuth. Using the disambiguation algorithm described, the next shots, that is the shots that have a failed azimuth less than 15 °, were disambiguously 95% of the time. The disambiguation algorithm produced correct results for more distant shots 75% of the time. The accuracy of disambiguation is assumed to vary depending on the geometry of the sensor array and the distribution

30 presumed of the shots, being easier to disambiguate the shots that have little elevation.

The solutions described for the trajectory of the projectile were obtained without rebuild detection. However, it has been found that even a weak cannon mouth signal, or a cannon mouth signal received only over a limited number of channels, can be advantageously used to improve the range of determination and disambiguation.

Figure 9 schematically shows a diagram of an arrival time model (TOA), which is described in greater detail in US Pat. No. 6 178 141. The TOA model can be used to estimate the trajectory of the projectile and the direction of the shooter in relation to the location of the sensor. The TOA model is based on a

P (P, P, P) of a shooter; elxyz

ballistic model that takes into account certain physical characteristics, relative to the flight path of the projectile, such as air density (which is related to temperature); the position

azimuth and elevation angles of the rifle's mouth; the initial velocity of the projectile (or equivalently the Mach number); and the speed of sound (which varies with the temperature / density of the air). With this ballistic model, it is possible to accurately calculate the moment when the shock wave and the rebuild reach a specific point in space.

P

R

S

As described in the diagram of Figure 9, the shooter is located at a point relative to an origin (0, 0, 0), the various sensors are located at points

. The vector from the shooter to the jth sensor is

(Px, Py, Pz) with

(Sxj, Syj, Szj) and the trajectory of the

j

TO

D

bullet is shown coming out of the shooter in the direction of,

the closest approach point (CPA) of the bullet to the jth sensor is ||

D

| =

| sin (1), and the path followed from the point at which the shock wave is radiated from the path to the jth sensor is S (the j index of the sensors has been omitted). The Mach angle of the bullet is 8 = sin-1 (1 / M), M =

V / c0. M the number of Mach of the projectile, V is the supersonic velocity of the projectile and c0 is the speed (dependent on pressure and temperature) of the sound. The 'failed angle' between the path and the jth sensor is

twenty 1. The trajectory is characterized by its azimuth angle measured counterclockwise from the x axis, in the x - y plane, and by its elevation angle measured upwards from the x - y plane. The equations that define the arrival time of the shock wave tj and the unit vector in the jth sensor are described in terms of these geometric quantities.

TO

TO

S

The arrival time is equal to the time 25 the sound is radiated to the jth sensor, plus the time it takes for the shock wave to travel the distance

it takes from the projectile to travel the distance to the point in the

V

S

from the point of radiation to the jth sensor,

c0

where t0 is a temporary reference (activation time) and c0 is the speed of sound. The angle of Mach 8 is also indicated in Figure 9.

30 It can be safely assumed that the velocity V of the projectile remains constant over the distance corresponding to the sensor separation, so that there is an insignificant loss of speed between the times in which the projectile radiates to the different sensors. However, it is known that over long distances the projectiles slow down due to air resistance. The air resistance can be expressed by a coefficient of aerodynamic drag Cb that depends on the shape of the bullet and the caliber of the bullet. A model

35 Ballistic mathematics derived from physical principles can predict the arrival time of a shock wave at any general point in space, based on a complete set of parameters that describe the projectile (for example, by its coefficient of aerodynamic drag Cb, its velocity initial and surrounding air density, which are known in advance).

The parameters required for an exact calculation are typically unknown in a realistic configuration, such

40 Like a battlefield However, the scope estimation can be significantly improved by an iterative process, shown in the form of a process flow diagram 200 in Figure 10, which takes into account the deceleration of the projectile velocity along the trajectory. Process 200 begins step 202, with the following assumptions:

c0 = speed of sound modified by external temperature / air pressure (~ 340 m / s) Cb = nominal efficiency of aerodynamic drag, averaged over the intended weapons

V0 = initial projectile velocity when fired, averaged over the planned weapons

M0 = V0 / c = initial Mach number of the projectile

In step 204, a first estimate of the distance of the shooter D0 is calculated, using the difference between arrival times (TDOA) measured Tms, and an arrival angle at between the sound of the muzzle and the shock wave, in the array of sensors, and assuming an initial velocity V0 and a constant number of Mach M0, according to the equation

With these assumptions, the velocity of the projectile at a distance from the position of the shooter P, can be calculated in step 206 from the equation

so that the time in which the projectile travels the distance, a, along the long trajectory remains, step 208,

The angle 8 is related to the Mach Ma number by the equation

where the Mach Ma number is initially set to M0. It should be noted that the instantaneous velocity of the bullet is set at the speed of sound (i.e., Ma = 1) if the calculated velocity of the bullet becomes less than the velocity of the

sound. The revised distance a = | A | in step 210 it is then like

twenty The angles a, 1 and 8 are related by the equation (a + 1 + 8) = 90 °. Then, the process 200 returns to step 206 by inserting the calculated value for the distance a, in the above equations for Ma and Ta, respectively providing an updated Mach number Ma and a displacement time of the updated Ta bullet, for the distance traveled to. The TDOA measured Tms and the updated values ​​calculated for Ta ya, are then used to successively update the D value for the shooter's reach:

This process is repeated until a maximum number of iterations has been reached, or the scope value D converges, as step 212 is determined.

Process 200 also verifies in step 214 if the value of the revised scope D = | D | for the distance between the handle and the sensor set, it is a "reasonable" value, in which case the process 200 ends in step 216. For example, the value for D can be considered valid if the distance a, traveled by the projectile , and the distance

between the sensor and the point at which the sound wave is radiated from the projectile to the sensor, they are

5 valid numbers, that is, it is not a NAN. A NAN is a special floating point value, which represents the result of a numerical operation that cannot return a valid numerical value, and is typically used to prevent the propagation of errors through a calculation. In addition, a and 1 must both be less than a predetermined threshold value, indicating that the projectile has undoubtedly been fired towards the sensor array.

As mentioned above, the pair of numbers (Tms, a) is initially used to calculate the scope of the

10 D0 shooter in the zero-th approach, disregarding changes in projectile speed along the trajectory. If the iterative process 200 described above does not return a consistent geometry that supports the pair of numbers (Tms, a), the solution is discarded.

Even if it is not possible to obtain an exact solution, one objective is to find values ​​for the reach of the D-handle and the failed azimuth and elevation angles (which are related to 1) that fit more closely to a measured TDOA wave of shock and a TDOA measured from the mouth of the barrel. As already mentioned, the TDOA shock wave only between the various sensors can be measured reliably in most situations. The azimuth of the shooter and the elevation of the shooter, but not the scope of the shooter, can be determined from the TDOA only of the shock wave, using the known coordinates of the sensor array (Sxj, Syj, Szj). It will be assumed that TDOA Tms can also be measured between the shock wave detected and the sound of the mouth

twenty of the barrel, so that the sound of the barrel's mouth may not be detected by all sensors.

If in step 214 it is determined that the iterative process 200 does not return a valid result, then the process 200 attempts to calculate the scope of the shooter, invoking an evolutionary genetic algorithm (GA) 300. The GA mimics the natural evolutionary principles, and applies them to search and optimization procedures. A GA begins its search with a random set of solutions, rather than just a solution. Once it has

25 created a random population of solutions, each one is evaluated in the context of the non-linear programming problem, and an aptitude (relative merit) is assigned to each solution. In one embodiment, the aptitude can be represented by the Euclidean distance between a calculated solution and the measured solution, for example by

Intuitively, an algorithm that produces a lower value of / Tmin is better.

30 A schematic flow chart of the GA 300 process is shown in Figure 11. The process 300 uses a difference between arrival times (TDOA) Tms and the arrival angle as previously measured by process 200, step 302. By way of example , you define a number of 3-tuples that have the values ​​{scope, MA, ME}, as

initial population at stage 304, where SCOPE is the reach of the shooter D = | D | shown in figure 9, MA is the failed azimuth and ME is the failed elevation. The MA and ME values ​​indicate to what extent the bullet has failed the

35 white, in azimuth and lifting space. In the illustrated example it is assumed that the target is the sensor matrix. In step 304, the initial population was created by random selection of the 3-tuple, covering a significant and reasonable range of values:

Shooter Range = {1000, ..., 3000} [meters],

Azimuth Failed = {-20, ..., 20} [degrees],

40 Failed Elevation = {-20, ..., 20} [degrees].

The calculation follows a process similar to the one outlined above for the shock wave solution only. Initially for the

Generation Gen = 0, step 306, the position vector of the shooter P (Px, Py, Pz) is calculated for each 3-tuple with the azimuth and elevation of the previously determined handle, and an assumed SCOPE for 3-tuple. Assuming

an initial number of Mach M0, the vector A (Ax, Ay, Az), that is to say the position from which the sound of the

wave, it is calculated as the MA and ME values ​​for each 3-tuple, step 308. The distances are also calculated

D = SJ -P between the shooter and each sensor j that detects a shock wave.

For each 3-tuple, angle 1 is calculated from the equation

where the symbol "." Indicates the scalar product of two vectors. The updated values ​​for distance a, the travel time Ta of projectile 5 over distance E, and the number of Mach Ma, are calculated by inserting the value calculated for 1, and the values ​​initially assumed for Ma = M0 already, in the above equations for Ma, Ta, a and D, step 312. This process is iterated several times for each of the 3-tuples, for example 3 times as indicated in step 312, after which the residue / Tmin defined above, which includes the cannon mouth signal, is calculated for every 3-tuple in the stage

314.

10 In step 316 it is verified whether a maximum number of iterations has been reached for the GA, for example 25 iterations. If the maximum number of iterations has been reached, then process 300 stops at step 320, returning the 3-tuple with the minimum residue. Otherwise, process 300 creates a new population through a crossover and mutation operation, step 318, and the generation counter is increased by one, step 322.

In each generation the "best" individual is allowed to survive without mutation, while the 100 individuals

fifteen Better estimating their aptitude also survive, but they are used to create the next 100 individuals from pairs of these survivors, with the crossing / mutation operators listed in the following table 2.

As an example, the following crossover and mutation operators were used to show process 300:

Table 2

<dl><dt>Name of the operator </dt><dd>Type of Operator Probability Description </dd></dl>

<dl><dt>Crossing - Failed Azimuth </dt><dd>Crossing 0.5 Swap failed azimuth between two chromosomes </dd></dl>

<dl><dt>Crossing - Failed Elevation </dt><dd>Crossing 0.5 Exchange failed elevation between two chromosomes </dd></dl>

<dl><dt>Crossing - Scope Failed </dt><dd>Crossing 0.5 Swap range between two chromosomes </dd></dl>

<dl><dt>Mutation-Field </dt><dd>Mutation 0.3 Replace a given field (with a probability of 0.25 per field) with a new randomly selected value, within the range </dd></dl>

<dl><dt>-Incremental mutation </dt><dd>Mutation 0.4 Introduce small disturbances in all fields of a chromosome (within ± 2 meters for the shooter's reach; within ± 0.1 ° for azimuth and failed elevation) </dd></dl>

twenty The GA 300 process is executed with an initial population of 200 different 3-tuples, with a fill rate of 50, for a total of 25 generations. The GA is executed 5 times in parallel with different sets of initial 3-tuples, and the solution with the lowest residue is selected, as the final solution for REACH, failed azimuth and

failed lifting of the handle, which allows the calculation of a vector D.

Recent experimental trials have indicated a decrease in ambiguous shots from 95% to 8% over the

25 same set of data, through the use of at least one signal channel from the muzzle, in addition to 5 or more shockwave channels, which is a significant improvement over solutions only with shockwave.

The calculations do not take into account that the deceleration of the projectile along its trajectory tends to overestimate the scope. For certain geometries and shots sufficiently distant, this overestimation may exceed 20%. The process described above eliminates this bias with respect to the estimation of the scope, for

30 long range shooting detections.

As described above, often ambiguous solutions only with a shock wave can be disambiguated by comparing the residues from two different paths, and selecting the trajectory with the least residue.

If the rebufo signals are detected on 4 or more sensor channels, then the wave algorithms of

5 The muzzle-cannon described above can be used to unambiguously determine the location of the shooter, regardless of the number of shockwave channels. If the rebufo signals are detected in less than 4 sensors, but the shock wave signals are detected in 5 or more shock wave channels, then the GA mentioned above can be used with a residue or modified cost function, by what any available cannon mouth signals are "mixed" in the optimization function, for

10 disambiguate the solution only with a shock wave and / or refine the estimate of the range of the shooter. However, if less than 3 cannon mouth channels and less than 3 shock wave channels are detected, then an alert can be activated without attempting to locate the shooter.

The signal from the muzzle may not be detected reliably on all channels, because:

1. The level of detection on one or more channels is too low for detection with confidence.

fifteen 2. The energy from the muzzle is not discernible in the raw signal, causing the system to correlate with 'noise', providing unreliable TDOA estimates.

3. The echoes from the shock wave may be larger than the rebuild, and may arrive before the rebuild causing the system to falsely detect the shock wave like the cannon mouth.

With a rebufo signal detected only on some channels, in this situation the residue can be defined as

where the first term for the rebufo is added on the reduced number of sensors (&lt;4) that detect the rebufo, and j is added on the sensors that detect the shock wave (typically, all the sensors).

As demonstrated by the examples described above, the rebufo signal provides information

25 important on the azimuth of the shooter, and therefore on the trajectory of the projectile, in comparison with a single shock wave solution, so that the solution of the calculated trajectory aligns better with one of the ambiguous solutions, that is, the solutions are therefore disambiguous.

Without at least some reliable barrel mouth signals, a significant number of ambiguous solutions can only be generated with a shock wave, especially at long distances from the shooter, which is less desirable than a

30 smaller number of solutions not ambiguous but less precise.

In the case of a potentially unattainable cannon mouth detection, initially the attempt can still be made to detect cannon mouth signals, for example to find a trace of a noisy signal, and calculate the resulting TDOA. The detection of the muzzle is estimated as reliable if signals from the muzzle are found on a sufficient number of sensors, with sufficient cross-correlation between the channels,

35 and if there is a sufficiently strong correlation between the signal of the cannon mouth and the corresponding raw band on each channel (a series of periods displaced, to take into account filtering delays).

In another case, at least the muzzle signals that have shown insufficient correlation are erased, and the following logic of 'coarse barrel detection' is invoked:

<dl><dt>-</dt><dd>Find maximum shock wave energy after a shock wave. Label these maximums as probable 40 'echoes of the shock wave', thereby excluding them from the overflows.</dd></dl>

<dl><dt>-</dt><dd>Determine a maximum time that would take the wave of the canyon to pass through the sensor array, and define a "window" that has the corresponding duration. Find maximum energy from the muzzle, by moving this window substantially through all channels of the detector after the detected shock wave, jumping sections in the detected signal that have been identified as echoes of the shock wave. Integrate the</dd></dl>

Four. Five energy over the window, that is, look for the maximum of:

It represents a measure of the energy, for example of the rebufo, measured by the nth sensor. The term (i + j) indicates the detection channel, with i indicating a discrete time interval between the moment the shock wave was detected and the beginning of the window, and j represents a time interval measured from the beginning of window.

where the square of

To discriminate against noise, the maximum energy in the window that produces the maximum of the function fmax (i) is verified, to determine if the peak of energy in the maximum is less than the energy through the window, in a given relationship factor. If this is the case, then the signal in the window is identified as a detection of the muzzle, and cross correlation is carried out on all the channels in the band of rebufo mb, to determine the TDOA of the muzzle .

Then, the determined rebufo signal can be used to determine the range of the shooter and / or to disambiguate the shock wave signal, as described above.

In summary, the described system can accurately, quickly and often unambiguously provide the direction of the shooter and the trajectory of the bullet, depending on only shock wave measurements. The disambiguation can be improved, and the range of the shooter can be estimated if at least one weak waveform of the rebufo is also detected. The system is relatively insensitive to false indications of the shooter in response to vehicle vibration and noise, such as wind noise, firecrackers or close shots in directions away from the system.

It should be mentioned that the system that detects the shock wave signals, performs two tests on the initial waveforms, to determine if in fact the signal can be attributed to the shock wave. First, the total energy measured in a frequency band between approximately 700 Hz and 10 kHz is compared with an empirical threshold value. Only if this threshold value is exceeded, can the signal form be considered as coming from a shock wave. Second, the temporal extension of the initial positive pressure peak, detected, must be greater than approximately 70 μs and less than approximately 300 μs. These criteria provide the system with immunity against impulsive noise, such as firecrackers and non-threatening shots. If these tests are not passed, the detected waveform is not considered a shock wave, and a solution for the shooter is not tested.

While the invention has been disclosed in relation to the preferred embodiments shown and described in detail, various modifications and improvements can be made thereon. By way of example, although illustrative embodiments are described having acoustic sensors such as microphones, this does not have to be the case. Instead, other types of mechanical or electrical pressure sensitive sensors can be used. In addition, the values ​​provided in Tables 1 and 2 for various operators are provided as examples only, and other values ​​may be selected depending on the actual field conditions.