EP2196973A1 - Traffic information unit, traffic information system, vehicle management system, vehicle, and method of controlling a vehicle - Google Patents

Abstract

A traffic information unit (MD1, MD2, MD3) according to the invention comprises a facility (MI) for tracking vehicle state information of individual vehicles present at a traffic infrastructure and a facility (T) for transmitting said vehicle state information to a vehicle (70B, 70E). A traffic information system may comprise a plurality of these traffic information units. The invention further comprises a vehicle management system (C) for a target vehicle (70B, 70E) that is capable of receiving and using the vehicle state information and a vehicle provided therewith is provided. Additionally a method for controlling traffic is provided.

Description

BACKGROUND OF THE INVENTION Field of the invention
• The present invention relates to a traffic information unit.
• The present invention further relates to a traffic information system.
• The present invention further relates to a vehicle management system.
• The present invention further relates to a vehicle provided with a vehicle management system.
• The present invention further relates to a method of controlling a vehicle.
Related Art
• Cruise control systems that maintain the speed of a target vehicle at a predetermined velocity are well-known. More recently adaptive cruise control systems were developed that also adapt the speed of the target vehicle to the state (e.g. relative position and speed) of a lead vehicle, directly in front of the target vehicle. The (partial) relative motion state of the lead vehicle is for example determined by radar measurements. Still more recently cooperative cruise control systems were developed that not only take into account the state of the lead vehicle but also from the state of one or more vehicles in front of the lead vehicle.
• Cooperative cruise control has the potential to improve traffic safety as well as traffic flow, as the control system can better anticipate the traffic situation than an adaptive cruise control system. In case of traffic of vehicles only provided with adaptive cruise control a sudden breaking of one of the vehicles tends to cause a shock-wave, as each vehicle only changes its state in response to the change of state of its immediate predecessor. Contrary thereto, a target vehicle provided with a cooperative cruise control system can also react to a change in state of another vehicle not directly leading the target vehicle provided with a cooperative cruise control system. This allows the target vehicle to more gradually adapt its state, e.g. its velocity. This is favorable for traffic flow and traffic safety. It is however a drawback of this system that it is dependent from input data from the motion state estimator mounted at other vehicles in the neighborhood (the motion state estimator is typically part of the CACC installed, thus this input is available via the cooperation between CACCs (hence the name)). Although the implementation of a cooperative cruise control system potentially allows for an improved traffic safety and traffic flow, such a control system is unreliable unless a relatively high fraction of the vehicles is provided with such a control system. Accordingly there is a need to provide a more reliable solution for improving traffic safety and traffic flow.
SUMMARY OF THE INVENTION
• According to an aspect of the invention a traffic information unit associated with a traffic infrastructure is provided comprising
• a facility for tracking vehicle state information of individual vehicles present at the traffic infrastructure,
• a facility for broadcasting said vehicles state information to other vehicles at the traffic infrastructure.
  A traffic information unit is considered associated with a traffic infrastructure if it has a sensor system using sensors that are mounted at an at least substantially fixed position related to the traffic infrastructure. For example the sensor system may comprise sensor nodes that are embedded in the traffic infrastructure. In order to allow a tuning of the sensor system the sensor nodes may arranged movably at a holder that has a fixed position with respect to the traffic infrastructure.
  In an embodiment the traffic information unit further comprises
• a sensor system comprising a plurality of sensor nodes for sensing vehicles arranged in the vicinity of a traffic infrastructure for carrying vehicles,
• communication means coupled to the sensor system, wherein the facility for tracking uses information communicated by the sensor system.
• According to a further aspect of the invention a vehicle management system is provided for target vehicles comprising a communication system arranged for receiving vehicle state information relating to surrounding vehicles from a traffic information unit, inputs for receiving state information from the target vehicle and a control system for providing control signals for controlling a state of the target vehicle using the other vehicles' state information retrieved from the traffic information system and the motion state of the target vehicle.
• According to a further aspect of the invention a vehicle with such a vehicle management system is provided.
• The traffic control system according to the present invention comprises a traffic information system that builds and maintains a real-time database of all vehicles currently using a traffic infrastructure. This enables a vehicle provided with a vehicle control system to receive status information of vehicles in its environment. In an embodiment said status information is only provided upon request. This allows for a power reduction as the transmitters do not have to be active when no such requests are received. Alternatively the transmitters may be active permanently and transmit this information unconditionally on a unidirectional basis. This is favorable if a large number of vehicles instrumented with a vehicle management system is present. In an embodiment the traffic information unit may have a first mode wherein vehicle status information is only transmitted upon request, e.g. when a low traffic density is detected and a second mode wherein the vehicle status information is permanently transmitted, e.g. during rush hours.
• Basically the traffic information system may broadcasts vehicle state information for the part of the infrastructure observed by the traffic information system. If desired the information may be restricted to information related to vehicles within a predetermined radius of a transmitter.
• Information to be transmitted may include not only vehicle state information relating to the lead vehicle (i.e. the vehicle directly in front of the target vehicle), but also vehicle state information relating to other vehicles in front of the lead vehicle that could not be observed by an on-board radar system. Also vehicle state information relating to vehicles behind the target vehicle may be included in the query set. As the traffic control system can provide status information, not only of the lead vehicle, but also of other vehicles in front of the target vehicle, the vehicle control system can better anticipate for events occurring at the road in front of the target vehicle, allowing for a smoother and safer control. It is not necessary that a large fraction of vehicle at the traffic infrastructure is provided with a vehicle control system according to the invention. The vehicle control system of each instrumented vehicle will operate reliably using the information transmitted by the traffic information system. Each of these instrumented vehicles can use the full vehicle map provided by the traffic information system according to the present invention and therewith reliable adapt its own motion to the If only a relatively modest fraction of the vehicles present at the road is provided with the inventive vehicle control system, these vehicles will already act as a buffer for smoothing traffic flow. The smoother traffic flow allows for a reduction in fuel consumption and air pollution. This would not be the case if the same number of vehicles were provided with a cooperative cruise control system, as the functioning of the cooperative cruise control system relies on the presence of other vehicles having the same cooperative cruise control system. Moreover as it is guaranteed that the necessary traffic data is provided by the traffic information system coupled to the traffic infrastructure, it becomes more attractive for owners of vehicles not provided with a vehicle management system according to the invention to achieve such a vehicle management system.
• Other applications of the present invention are possible. One of them is formation driving. Because the traffic information system provides the vehicles instrumented with a vehicle control system with state information in its environment, and therewith allows the vehicle control system to anticipate for events ahead of the target vehicle, the vehicle control system can maintain short distances to its predecessor.
• Incident management is a further example. The traffic management system can provide information to a target vehicle about incidents ahead of the target vehicle and enforce safety measures. The safety measures may include a gradual braking of the target vehicle, a deviation of the target vehicle to an alternative route, a warning to the driver of the vehicle and/or a warning to other drivers by light signals.
• In an embodiment of the traffic information unit the sensor system comprises a plurality of sensor nodes that each provides a message indicative for an occupancy status of a detection area of a traffic infrastructure monitored by said sensor node. The traffic information system further comprises at least one message interpreter that includes:
• a vehicle database facility comprising motion state information of vehicles present at the traffic infrastructure, the state information of the vehicles including at least the vehicle position,
• an association facility for associating the messages provided by the sensor nodes with the state information present in the vehicle data base facility,
• a state updating facility for updating the state information on the basis of the messages associated therewith.
• In the traffic information unit according to this embodiment vehicles can be tracked with relatively simple and cheap means. It is sufficient that the sensor nodes merely provide a message that indicates whether a detection area associated with the sensor node is occupied by a vehicle or not. This makes it economically feasible to apply the traffic information unit to large traffic infrastructures.
• Suitable sensor elements for use in a sensor node are for example magnetic loop sensors, magneto restrictive sensors. These sensor elements determine whether their associated detection area is occupied by detection of a perturbation of the earth magnetic field.
• Preferably each sensor node is provided with a wireless transmission facility that transmits the sensed data to a receiver facility coupled to the association facility. This facilitates installation of the sensor nodes. Furthermore it is attractive if the sensor nodes provide their message at an event basis, e.g. if a perturbation of the earth magnetic message interpreter and minimizes power consumption of the sensor nodes.
• In an embodiment a sensor node may have a set of sensor elements that are clustered in a sensor module. Such a sensor module is for example a camera that monitors a part of the traffic infrastructure, wherein each photosensitive element of the camera serves as a sensor element of the vehicle tracking system. A camera may be used for example if a perturbation of the earth magnetic field can not be measured. This is the case for example if (parts of) the infra structure comprises metal components e.g. a bridge.
• It is not necessary that the detection areas of the sensor elements are complementary. The detection areas may overlap, or spaces may exist between the detection areas. It is sufficient that the detection areas have a scale that is smaller than the vehicle to be tracked, e.g. a size of at most 1 m² and a maximum diameter of not more than 1 m.
• Surprisingly it has been found that it is advantageous if the sensor elements are randomly distributed over the traffic infrastructure. As compared to an arrangement wherein the sensor elements are regularly distributed with the same average number of sensor elements per unit of area, a more accurate estimation of the state of the vehicles can be obtained.
• Independent traffic information units are particularly suitable for providing vehicle state information for relatively small traffic infrastructures. Particularly for larger traffic infrastructures a traffic information system is provided that comprises at least a first and a second traffic information unit according to the present invention. The first and the second Traffic information unit are associated with mutually neighboring sections of the traffic infrastructure and are arranged to mutually exchange state information.
• In this way a traffic information system is provided that can be easily expanded with one or more additional traffic information units if required. A new traffic information unit needs only to communicate with the traffic information units arranged for neighboring sections. For example if a certain road is already provided with a traffic information system, it is sufficient to provide for a communication facility between the information unit for the last section of said traffic information system and the new traffic information unit for the appended section. As the traffic information units, merely exchange state information and not the unprocessed messages from the sensor nodes the amount of communication between the traffic information units is modest.
• An embodiment of a vehicle management system further comprises communication means for exchanging vehicle state information with surrounding vehicles and a selection facility for selecting one or more of vehicle state information obtained from the surrounding vehicles and information received from the traffic information system as the vehicle state information to be used by the control system.
• The selection made by the selection facility may for example depend on the availability of reliable information. For example in an area where the traffic infrastructure is provided with a traffic information system, the selection facility may automatically select the traffic information system as the source of state information. In an area where no traffic information system is present, it may select the information provided by surrounding vehicles. Alternatively the selection may be more fine grained. It may select for example to receive velocity information from the surrounding vehicles themselves if such information is available and to receive the remaining information from the traffic information system.
• According to a further aspect of the invention, a method of controlling traffic is provided, comprising the steps of
• observing vehicles from a fixed position,
• communicating the observations,
• tracking motion states of individual vehicles using the communicated observations
• transmitting said information about said tracked states to a vehicle instrumented with a vehicle management system.
BRIEF DESCRIPTION OF THE DRAWINGS
• These and other aspects are described in more detail with reference to the drawings. Therein:
• Figure 1 schematically shows a spatial arrangement for various components of a traffic information system according to the invention,
• Figure 2 shows a functional relationship between various units of the system of Figure 1,
• Figure 3 shows another schematic view of the traffic information system,
• Figure 4 schematically shows a part of the traffic infrastructure provided with a plurality of sensor elements,
• Figure 5 schematically shows an overview of interactions between a traffic information system for an traffic infrastructure and a vehicle management system of vehicles using the traffic infrastructure,
• Figure 6 shows an embodiment of a vehicle management system according to the invention,
• Figure 7 shows an example of an embodiment of a sensor node in a traffic information system according to the invention,
• Figure 8 shows a possible hardware implementation of the sensor node of Figure 7,
• Figure 9 schematically shows a method performed by the sensor node of Figures 7 and 8,
• Figure 10 shows a message interpreter in a traffic information system according to the invention,
• Figure 11 shows a possible hardware implementation of the message interpreter of Figure 10,
• Figure 12 schematically shows a part of a traffic infrastructure,
• Figure 13 shows an overview of a method performed by the message interpreter,
• Figure 14 shows a first aspect of the method in more detail,
• Figure 15 shows a second aspect of the method in more detail.
• Figure 16 shows an example of an object to be detected at a reference position and orientation and at a different position and orientation,
• Figure 17 shows a definition of a set S and the equidistant sampled set Λ,
• Figure 18 shows detection of an object at multiple detection points,
• Figure 19 shows a definition of the set O_n of possible positions oⁱ _k for a single detection point,
• Figure 20 shows a definition of the set ON of possible positions oⁱ _k for multiple detection points,
• Figure 21 shows a derivation of O_N(q) given 2 detections and 2 different samples of q,
• Figure 22 shows (Left) determination of -Λ, (right) the object's possible position set ^On given d_n and q,
• Figure 23 shows (left) the mean of all Gaussians from f(o| z₁, θ) and f(o|z₂, θ); (right) The selection of means of the Gaussians from f(o|z₁, θ) and f(o|z₂, θ), of which their mean ô_ni is close or in C_N(θ),
• Figure 24 shows an association result with event-based data-association,
• Figure 25 shows an association result with Nearest Neighbor data-association,
• Figure 26 shows time sampling of a signal y(t),
• Figure 27 shows event sampling of a signal y(t),
• Figure 28 shows event sampling: Send-on-Delta,
• Figure 29 shows the Gaussian function,
• Figure 30 shows a top view of the Gaussian function,
• Figure 31 shows an approximation of Λ _Hke (y_n ) as a sum of Gaussian functions,
• Figure 32 shows position, speed and acceleration of a simulated object,
• Figure 33 shows a position estimation error for various methods,
• Figure 34 shows a speed estimation speed for various methods,
• Figure 35 shows a factor of increase in estimation error after z_ke , or y _ka .
DETAILED DESCRIPTION OF EMBODIMENTS
• Figure 1 and 2 schematically show a traffic information system comprising a plurality of traffic information units. Therein Figure 1 schematically shows how in an example embodiment various components of the system are arranged. Figure 2 shows a functional relationship between various units of the system. As shown in Figure 1, the traffic information units comprise a sensor system with a plurality of sensors (indicated as black dots) for sensing vehicles (indicated by open hexagons) arranged in the vicinity of a traffic infrastructure 80 for carrying vehicles. The sensors are provided with communication means to transmit sensed information to a facility MI for identifying and tracking states of individual vehicles using information communicated by the sensors. Although in this embodiment the sensors are only capable of transmitting information towards the facilities MI, in another embodiment, they may also be capable of bidirectional communication. In that embodiment sensors can form a network, that can guide the information in an indirect way to the facilities MI. In this embodiment each of the facilities MI is responsible for monitoring a particular section 80A, 80B, 80C, 80D of the infrastructure 80. In Figure 1 and 2 only four facilities MI are shown for clarity.
• Figure 3 is another schematic view of the traffic information system. Figure 3 shows how sensor nodes 10 transmit (detection) messages D to a message interpreter MI in their neighborhood. The message interpreters MI may also communicate to each other via a communication channel 60 to indicate that a vehicle crosses a boundary between their respective sections and to exchange a status of such a vehicle. As shown in Figure 3, the traffic information system comprises a plurality of traffic information units MD1, MD2, MD3. Each traffic information unit MD1, MD2, MD3 comprises a respective subset of the plurality of sensor nodes 10 for monitoring a respective section of the traffic infrastructure and a respective message interpreter MI. The traffic information system further has a communication facility 60 for enabling traffic information units MD1, MD2, MD3 of mutually neighboring sections to exchange state information. The traffic information system further comprises client information modules CIM for providing status information related to the infrastructure 80. The status comprises for example statistical information, such as an occupation density and an average speed as a function of a position at the traffic infrastructure 80. The facilities MI and the client information modules CIM are coupled to each other via a communication backbone. This allows the client information modules CIM to request said information for arbitrary regions (indicated by dashed boxes) of the infrastructure 80 that may extend beyond the boundaries for individual facilities MI.
• Figure 4 schematically shows a part of the traffic infrastructure that is provided with a plurality of sensor nodes j having position c_j. The sensor nodes have a detection area with radius R. A vehicle i is present at the infrastructure having a position (vⁱ _x, vⁱ _y). In this case if the vehicle substantially covers the detection area, e.g. more than 50%, the sensor node sends a message D that the detection area is occupied (indicated in gray). Otherwise the sensor node sends a message that the detection area is not occupied (white).
• As shown in Figures 3 and 5, the traffic information system is further provided with a facility T for transmitting state information derived by the traffic information system to a particular vehicle upon request. Each transmitter T has a transmission range TR. Preferably the transmission ranges of the transmitters together define a continuous area having a substantial length and over a full width of the infrastructure where state information is available. A plurality of transmitters may be coupled to each traffic information unit MD1, MD2, MD3. Preferably the transmitters T selectively transmit vehicle state information related to vehicles within their transmission range and optionally in a neighborhood thereof.
• As shown in Figure 5 some 70B, 70E of the vehicles 70A,...,70E present at the traffic infrastructure 80 are provided with a vehicle management system C. The vehicle management system C comprises a communication system R arranged for receiving vehicle state information relating to surrounding vehicles from the traffic information system, e.g. here from the traffic information unit MD1. The traffic information unit MD1 transmits the motion state of the surrounding vehicles to the target vehicle (e.g. 70B) provided with a vehicle management system C, using the wireless link between the transmitter T and the communication system R of the vehicle management system C. This information is stored in a local vehicle status data base C0. The vehicle management system C further has inputs C1 for receiving state information from the target vehicle 70B. The state information may include information related to a momentaneous position, e.g. obtained by GPS, speed obtained by GPS or using odometry, an acceleration derived by odometry or by an inertial sensor and a direction e.g by using a compass or a by a gyro.
• A control system C2 uses this information in the local vehicle status database C0 and the state information received at inputs C1 to provide control signals at output C3 for controlling a state of the target vehicle, e.g. a speed or an orientation of the target vehicle (70B).
• In the embodiment shown, the vehicle management system C also has an bidirectional link C4 for additional communication purposes. This link can be used to negotiate and coordinate actions among vehicles (e.g. requesting/granting free space, joining/leaving platoon, etc.). The system C further has an input C5 for receiving user control commands. This allows the user to set an authorization level, i.e. control the extent to which the system C controls the vehicle, e.g. the user may allow the system only to provide warnings, may allow the system to regulate a speed, to break the vehicle up to a predetermined maximum deceleration, and to control a travelling direction. In the latter case a user may for example instruct the system to carry out certain maneuvers, e.g. a merging between a sequence of vehicles in a neighboring lane.
• In the embodiment shown a further input C6 is present to receive navigation information. This information may be used for global control. For example dependent on a particular route to follow as indicated by the navigation information, the control system C2 may control the vehicle to another lane, taking into account the state of neighboring vehicles in local vehicle status data base C0.
• Output C7 may provide the user information about the current authorization level, about a current activity of the system C, to show warnings, and to request for input. The C7 output represents a man-machine interface and may be implemented in any form; it may use auditory, visual or sensory channels.
• An embodiment, wherein the traffic information system only provides the state information of neighboring vehicles upon request has the advantage that power is saved during intervals that no information is requested. Alternatively however, the transmitters T may permanently transmit the information relating to the vehicles present in its neighborhood.
• As the traffic information system can provide status information to an instrumented vehicle, e.g. 70B, not only of the lead vehicle 70C, but also of other vehicles 70D, ..., 70E in front of the target vehicle 70B, the vehicle control system C can better anticipate for events occurring at the road in front of the target vehicle 70B. This allows for a smoother and safer control. For example the traffic information system will also transmit the status information of vehicle 70D, indicating that this vehicle intends to change from the rightmost lane to the middle lane of the traffic infrastructure 80. Using this information, the traffic control system C of vehicle 70B may respond more gradually to the maneuver of vehicle 70D, than would be the case if vehicle 70B had only a simple cruise control system that merely responds to the behavior of a vehicle immediately in front. It is not necessary that a large fraction of the vehicle present at the traffic infrastructure is provided with a vehicle control system according to the invention. The vehicle control system of each vehicle will operate reliably using the information transmitted by the traffic information system. If only a relatively modest fraction of the vehicles present at the road is provided with the inventive vehicle control system, these vehicles will already act as a buffer for smoothing traffic flow. This can be illustrated by way of the following example. Presume that the vehicles 70A, ..., 70E are driving in the same lane, and that none of the vehicles 70A, ..., 70E is instrumented with a vehicle control system or is only instrumented with an adaptive cruise control system. In that case a sudden breaking of vehicle 70E would result in a shock effect that ripples through the chain of vehicles. However, even if only a part of the vehicles is instrumented with a vehicle management system according to the present invention say 70B, the situation is different. In that case, substantially at the moment that vehicle 70E breaks, this change in vehicle status information 70E is detected by the traffic information system and communicated to the vehicle 70B instrumented with vehicle management system according to the present invention. This allows vehicle 70B to anticipate for the shockwave that ripples through the sequence of vehicles 70C, 70D, 70E preceding it. Therewith the control system C2 of vehicle 70B can initiate a smooth breaking procedure starting substantially at the moment of the sudden breaking of vehicle 70E. This not only has positive consequences for the vehicle 70B itself, but also for the vehicles 70A behind it.
• The set of vehicles for which vehicle status information is transmitted by a transmitter T in the neighborhood of a target vehicle, e.g. 70B may include vehicles 70C,...,70E, may additionally or alternatively include vehicles 70A behind the target vehicle 70B. This vehicle status information may be used by the control system C2 to of vehicle 70B to moderate a breaking power of said vehicle 70B to prevent that a collision occurs with a vehicle 70A succeeding it.
• Figure 6 shows a further embodiment of a vehicle management system C according to the invention. Parts therein corresponding to those in Figure 5 have the same reference. The vehicle management system of Figure 6 further comprises communication means R1 for exchanging vehicle state information VS2 with surrounding vehicles. The vehicle management system C shown therein further comprises a selection facility SL for selecting one or more of vehicle state information VS2 obtained from the surrounding vehicles and vehicle state information VS1 received from the traffic information system as the vehicle state information VS to be used by the control system C2. The control system C2 further receives state information from the target vehicle (ST).
• The selection made by the selection facility SL may for example depend on the availability of reliable information. For example in an area where the traffic infrastructure is provided with a traffic information system, the selection facility may automatically select the state information VS1 provided by said traffic information system as the source of state information VS. In an area where no traffic information system is present, it may select the information VS2 provided by surrounding vehicles. Alternatively the selection may be more fine grained. For example it may select for example to receive velocity information from the surrounding vehicles themselves if such information is available and to receive the remaining information from the traffic information system.
• In the sequel an embodiment of a traffic information system is described. Therein Figure 7 shows an example of a sensor node 10. The sensor node 10, shown in Figure 7, is an assembly of a sensor element 12, a processing unit 14 (with memory), a clock-module 18 and a radio link 16.
• The sensor element 12 is capable of sensing the proximity of the vehicles to be tracked. The processing unit 14 determines if an object (vehicle) is present or absent on the basis of the signals from the sensor element 12. If an occupancy status of the detection area of the sensor changes, the processing unit 14 initiates a transmission of a message D indicating the new occupancy status and including a time stamp indicative of the time t at which the new occupancy status occurred. The message D sent should reach at least one message interpreter MI. A concrete implementation of the sensor node 10 is used for road vehicle tracking: in this case the sensor element 12 is a magnetoresistive component, which measures the disturbance on the earth magnetic field induced by the vehicles. Alternatively, a magnetic rod or loop antenna may be used for this purpose.
• Figure 8 shows a possible implementation of the hardware involved for the sensor node 10 of Figure 7. The sensor element 12 is coupled via an A/D converter 13 to a microcontroller 14 that has access to a memory 15, and that further controls a radio transmitter 16 coupled to an antenna 17.
• Figure 9 schematically shows a method performed by a sensor node to generate a message indicative for occupancy status of a detection area of the sensor node.
• Starting (Step S1) from an off-state of the sensor node, input from the A/D converter is received (Step S2). In a next step S3, offset is removed from the sensed value.
• In step S4 it is determined whether the occupancy state of the detection area as reported by the last message transmitted by the sensor node was ON (selection YES) (vehicle present in the detection range) or OFF (selection NO) (no vehicle present in the detection range. This occupancy state is internally stored in the sensor node.
• In the first case, program flow continues with step S5. In the second case processing flow continues with step S9. In step S5 it is determined whether a signal value v obtained from the A/D converter, and indicative for an occupied status of the detection area is below a first predetermined value T_L. If this is not the case program flow continues with step S2. If however the value is lower than said first predetermined value then program flow continues with step S6. In step S6 it is verified whether the signal value v remains below the first predetermined value T_L for a first predetermined time period. During step S6 the retrieval of input from the A/D convertor is continued. If the signal value v returns to a value higher then said predetermined value T_L before the end of said predetermined time-period then processing flow continues with step S2. Otherwise the value for the occupancy state is internally saved as unoccupied in step S7, and a message signaling this is transmitted in step S8.
• In step S9 it is determined whether the signal value v obtained from the A/D converter, and indicative for an occupied status of the detection area exceeds a second predetermined value T_H. The second predetermined value T_H may be higher than the first predetermined value T_L. If the signal value does not exceed the second predetermined value T_H program flow continues with step S2. If however the value is higher than said second predetermined value T_H then program flow continues with step S10. In step S10 it is verified whether the signal value v remains above the second predetermined value T_H for a second predetermined time period, which may be equal to the first predetermined time period. During step S10 the retrieval of input from the A/D convertor is continued. If the signal value v returns to a value lower then said predetermined value T_H before the end of said predetermined time-period then processing flow continues with step S2. Otherwise the value for the occupancy state is internally saved as occupied in step S11, and a message signaling this is transmitted in step S12.
• A message interpreter, shown in Figure 10 and 11, consists of a radio receiver 20, coupled to antenna 22, a processing unit 24 (with memory 28) and a network interface 65, as well as a real-time clock 26. The network interface 65 couples the message interpreter MI via the communication channel 60 to other message interpreters.
• As shown in more detail in Figure 10, the radio receiver 20 receives the binary "object present" signals D (with timestamp) from the sensor nodes 10 via the radio link and runs a model based state estimator algorithm to calculate the motion states of the objects individually (i.e. each real world object is represented in the message interpreter). The accuracy and the uncertainty of the estimation depends on the sensor density. For accurate object tracking it is preferred to have coverage of multiple sensors per object.
• The message interpreter MI has a vehicle database facility 32, 34 that comprises state information of vehicles present at the traffic infrastructure.
• The message interpreter MI further has a sensor map 45 indicative for the spatial location of the sensor nodes 10. Alternatively, the sensor nodes may transmit their location or their position could even be derived by a triangulation method.
• The message interpreter MI further has an association facility 40 for associating the messages D provided by the sensor nodes 10 with the state information present in the vehicle data base facility 32, 34. The association facility 40 may associate the messages received with state information for example with one of the methods Gating, Nearest Neighbor (NN), (Joint) Probabilistic Data Association ((J)DPA), Multiple Hypothesis Tracker (MHT) and the MCMCDA.
• The message interpreter further has a state updating facility 50 for updating the state information on the basis of the messages D associated therewith by the association facility 40. Once the messages D are associated with a particular vehicle the state of that vehicle in a local vehicle data base is updated by the update facility 50.
• In the embodiment shown a global map builder 65 may exchange this updated information with global map builders of neighboring message interpreters via network interface 60 (wired or wireless), for example to exchange the motion state of crossing objects.
• In the embodiment shown the microcontroller 24 of Figure 11 processes the received messages D. The memory 28 stores the local and global vehicle map and the sensor map as well as the software for carrying out the data estimation and state estimation tasks. In an alternative embodiment separate memories may be present for storing each of these maps and for the software. Likewise dedicated hardware may be present to perform one or more of these tasks.
• There is no communication or any other interaction between the objects tracked. The result of the processing (i.e. the estimation of the motion states of all sensed objects) is present in the memory of the message interpreters in a distributed way.
  Message interpreters may run additional (cooperative) algorithms to deduct higher level motion characteristics and/or estimate additional object characteristics (e.g. geometry).
• For applications in relative small area, e.g. a parking place, or a traffic node, the vehicle tracking system may comprise only a single traffic information unit. In that case the global map builder is superfluous, and local vehicle map is identical to the global vehicle map.
• In the embodiment shown in Figure 3, each message interpreter MI for a respective traffic information unit MD1, MD2, MD3 comprises hardware as described with reference to Figure 10 and 11.
• Operation of the message interpreter is further illustrated with respect to Figures 12-15.
• Figure 12 schematically shows a part of a traffic infrastructure 80 having sections R_j-1, R_j, R_j+1. By way of example it is presumed that a vehicle moves in a direction indicated by arrow X from R_j-1, via R_j, to R_j+1.
• Figure 13 shows an overview of a method for detecting the vehicle performed by the message interpreter for section R_j, using the messages obtained from the sensor nodes.
• In step S20 the method waits for a message D from a sensor node. At the moment that a message D is received, program flow continues with step S21, where the time t associated with the message is registered. The registered time t associated with the message may be a time-stamp embedded in the message or a time read from an internal clock of the message interpreter.
• In embodiments wherein messages are indirectly transmitted to a message interpreter, e.g. by a network formed by sensor nodes it is advantageous if the embeds the time stamp in the message, so that it is guaranteed that the registered time corresponds to the observed occupancy status regardless any delays in the transmission of the message.
• In step S22, it is verified whether the detection is made by a sensor node in a location of section Rj that neighbors one of the neighboring sections R_j-1 or R_j+1. If that is the case, then in step S23 the event is communicated via the communication network interface to the message interpreter for that neighboring section. In step S24 it is determined which vehicle O in the vehicle data base facility is responsible for the detected event. An embodiment of a method used to carry out step S24 is described in more detail in Figure 14. After the responsible object O is identified in step S25, i.e. an association is made with existing object state information, it is determined in Step 26 whether it is present in the section Rj. If that is the case, control flow continues with Step S27. Otherwise control flow returns to step S20, where the state of object O is estimated. A procedure for estimating the state is described in more detail with reference to Figure 15. In step S28 it is determined whether the state information implies that the vehicle O has a position in a neighboring region R_j-1 or R_j+1. In that case the updated state information is transmitted in step S29 to the message interpreter for the neighboring region and control flow returns to step S20. Otherwise the control flow returns immediately to Step S20.
• A method to associate a message D at time t, with an object O is now described in more detail with reference to Figure 14.
• In a first step S40, a vehicle index i is initialized (e.g. i=1). In a next step S41, the current state known for the vehicle with that index i is retrieved from the vehicle database facility. In the next step S42 a probability is determined that the vehicle O caused the detection reported by the message D at time t. The vehicle index i is incremented in step S43 and if it is determined in step S44 that i is less than the number of vehicles, the steps S41 to S43 are repeated. Otherwise in step S45 it is determined which vehicle caused the detection reported by the message D at time t with the highest probability. In step S46 the index of that vehicle is returned as the result if the method.
• A method to estimate (update the present estimation of) the state of a vehicle is now described in more detail with reference to Figure 15.
• In step S60 the messages D₁,...,D_n associated with vehicle O are selected?
• In step S61 a probability density function is constructed on the basis of the associated messages.
• In step S62 the current state So and time to for object O are determined.
• In step S63 it is determined whether the time for which the state S of the vehicle O has to be determined is less than the time t₀ associated with the current state S₀.
• If that is the case, the state S determined by the estimation method is the state update of S0 to t, performed in step S65. If that is not the case, the state S determined by the estimation method is the state update of S0 to S0 in step S64. What does it mean?
• It is noted that other methods are possible to track vehicle state information of individual vehicles. For example vehicles could be provided with a transponder that signals their momentaneous position to the traffic information system.
• In the claims the word "comprising" does not exclude other elements or steps, and the indefinite article "a" or "an" does not exclude a plurality. A single component or other unit may fulfill the functions of several items recited in the claims. The mere fact that certain measures are recited in mutually different claims does not indicate that a combination of these measures cannot be used to advantage. Any reference signs in the claims should not be construed as limiting the scope. Further, unless expressly stated to the contrary, "or" refers to an inclusive or and not to an exclusive or. For example, a condition A or B is satisfied by any one of the following: A is true (or present) and B is false (or not present), A is false (or not present) and B is true (or present), and both A and B are true (or present).
• More details relevant for the present invention are described in the following Annexes:
• A1: Estimation and association for multiple target tracking based on spatially, distributed detections
• A2: On Event Based State Estimation
A1: Estimation and association for multiple target tracking based on spatially, distributed detections
• Summary. In this paper we consider the multiple object tracking problem with event-based observations. For that we predefine a number detection points which are spatially distributed along the road. Whenever the edge of an object crosses one of the detection points, the position of that detection point together with the time of the event are received by our tracking algorithm. We assume that objects can cover multiple detection points and propose a method to estimate the object's position and orientation from these detections using the shape of the object. Beside that another method is designed which associates newly received detections with a known object. The objects are tracked with an event-based state-estimator that uses the estimated position and orientation, although its design is out of the scope of this paper. Finally our tracking algorithm is critically assessed in a simulation of a parking lot.
1 Introduction
• In multiple target tracking [1-3] one aims to track all the objects/targets, which are moving in a certain area. Three basic problems arise from tracking objects. The first one is how to measure the object's position. The second one is to associate a certain measurement with its correct object and the third one is a state-estimator to keep track of all the objects. This paper considers the first 2 issues when objects are not measured but detected.
• Consider a system in which objects are detected when they cross a predefined 'detection' point. These detectors are triggered by the event that the object's edge crosses its position. However, they cannot distinguish between the objects. This paper describes a method in which a new detection is associated with the object that most probable generated it. Also, a method is described which estimates the position and orientation of the object given the observations in position and time due to the detections. Other examples in which sensor-data is generated due to an event can be found in [4,5].
• This paper is organized as follows. Section 2 defines background knowledge such as the notation of (object) variables and functions that are used throughout this paper. After that the problem is formulated in section 3 together with existing methods. Section 4 describes the approach which is taken in the design. A more detailed description of the estimation and associated is presented in Section 5 and 6 respectively. Finally both methods are tested in a small application example presented in Section 6 and conclusions are drawn in section 7. But let's start with the background information.
2 Background
• In order to be clear about notations and variables this section describes those that can be found throughout this paper.
2.1 Variables
• 

  defines the set of real numbers whereas the set

  

  defines the non-negative real numbers.

  

  defines the set spanned by the vectors e_x and e_y , the point p:=x·e_x +y·e_y is shortly denoted as p = (x,y) ^T. The set

  

  defines the integer values and

  

  defines the set of non-negative integer numbers. The variable 0 is used either as null, the null-vector or the null-matrix. Its size will become clear from the context.
• Vector x(t) ∈

  

  is defined as a vector depending on time t and is sampled using some sampling method. The time t at sampling instant k ∈

  

  is defined as t_k ∈

  

  . The variables τ _k ∈

  

  , x_k ∈

  

  and x _0:k ∈

  

  are defined as: $${\tau }_k:=t_k-t_{k-1},$$

  

  $$x_k:=x\left(t_k\right),$$

  

  $$x_{0:k}:=\left(\mathrm{<figure-callout\; id="0"\; label="x\; t"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">x</figure-callout>}\left(t_{}\right)\left({}_0\right)x\left({\mathrm{<figure-callout\; id="1"\; label="t"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">t</figure-callout>}}_1\right)\cdots x\left(t_k\right)\right)\mathrm{.}$$

  
• The matrix A(t ₂ - t ₁) ∈

  

  depends on the difference between two time instants t ₂ > t ₁ and is shortly denotes as A _{t 2 -r 1}
2.2 Functions
• The transpose, inverse and determinant of a matrix A ∈

  

  are denoted as A^T, A ^-1 and |A| respectively.
• Let us define the probability of the random vector x ∈

  

  as the scalar Pr(x) ∈ {0,1} and the conditional probability of x given the vector u ∈

  

  as the scalar Pr(x|u) ∈ {0,1}. The probability density function (PDF), as defined in [6] section B2, of the vector x ∈

  

  is denoted as p(x) and the conditional PDF of x given u ∈

  

  is denoted as p(x|u). The expectation and covariance of x are denoted as E[x] and cov(x) respectively. The conditional expectation of x given a vector u is denoted as E[x[u]. The definitions of E[x], E[x|u] and cov(x) can be found in [6] sections B4 and B7.
• The Gaussian function, shortly noted as Gaussian, depending on vectors x ∈

  

  and u ∈

  

  and on matrix P ∈

  

  is defined as: $$G\left(xuP\right):{ℝ}^n\times {ℝ}^n\times {ℝ}^{n\times n}\to ℝ,=\frac{1}{{\left(2\pi \right)}^{n/2}\left|P\right|}{e}^{-0.5{\left(x-u\right)}^T{P}^{-1}\left(x-u\right)}\mathrm{.}$$

  
• If p(x) = G(x,u,P), then by definition it holds that E[x] = and cov(x) = P.
• Assume we have the set

  

  ⊂

  

  and the vectors x ∈

  

  and y ∈

  

  . Then the function ∥x - y∥ ∈

  

  is defined as the distance between vectors x and y. The function |〈x -

  

  〉| ∈

  

  is defined as the shortest distance from vector x to set

  

  : $$\begin{array}{ll}|\langle x-ℂ\rangle |& :{ℝ}^q\to ℝ,\\ & :=\mathit{min}\left(\Vert x-c\Vert \right),\forall c\in ℂ\mathrm{.}\end{array}$$

  
2.3 Object variables
• Assume there exist an object which is moving in the 3D world space. This object is observed with, for example, a camera or sensors in the road. Meaning that the object is projected to a 2D space, i.e.

  

  . If we assume that the shape of the projected object is constant and known, then we can draw a smallest, rectangular box around the object.
• For the box we define a position-vector o = (x,y) ^T ∈

  

  , equal to the center of the box, and an orientation-vector θ ∈

  

  . In the case of o = 0 and θ = 0 the corners of this box, as shown in Figure 16, are defined in the set

  

  : $${ℂ}_0:\left\{{\mathrm{<figure-callout\; id="1"\; label="c"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">c</figure-callout>}}_1{c}_2{c}_3{\mathrm{<figure-callout\; id="4"\; label="c"\; filenames="imgf0008.png,imgf0024.png"\; state="\{\{state\}\}">c</figure-callout>}}_4\right\},\mathrm{with}{c}_i\in {ℝ}_{\mathit{xy}}\mathrm{.}$$

  
• Notice that for an object having a certain o and θ the new corner-positions of the object's box are calculated with

  

  . For that a rotation matrix T ∈

  

  is used as defined in (5). An example of the object's box for a certain o and θ is graphically depicted in Figure 16. $$T_k=\left(\begin{array}{cc}\mathit{cos}\left({\theta }_k\right)& \mathit{sin}\left({\theta }_k\right)\\ -\mathit{sin}\left({\theta }_k\right)& \mathit{cos}\left({\theta }_k\right)\end{array}\right)\mathrm{.}$$

  
• Beside the positions o and θ each object also has a certain shape or geometry which covers a certain set of positions in

  

  , i.e. the grey area of Figure 16. This closed set is denoted with S ⊂

  

  and is defined as the union of the open set of the object's body S_B ⊂

  

  and the closed set of the object's edge S_E ⊂

  

  , i.e. S := S_B ∪ S_E. The set S is approximated by a set of sampled position-vectors Λ = [λ₁,λ₂,···, λ _K ], with λ _i ∈

  

  . To define the vectors λ _i we equidistant sample the rectangular box defined by

  

  using a grid with a distance r. Each λ _i is a grid point within the set S as graphically depicted in Figure 17.
• The aim is to estimate position, speed and rotation of the object in the case that its acceleration and rotational speed are unknown. Therefore the object's state-vector s(t) ∈

  

  and process-noise w(t) ∈

  

  are defined as: $$s\left(t\right):={\left(o\left(t\right)\theta \left(t\right)\frac{\mathit{\delta o}\left(t\right)}{\mathit{\delta t}}\right)}^T,w\left(t\right):={\left(\frac{{\mathit{\delta o}}^2\left(\mathrm{<figure-callout\; id="2"\; label="t\; \delta t"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">t</figure-callout>},\right)}{{\mathit{\delta t}}^{}}\frac{\mathit{\delta \theta }\left(t\right)}{\mathit{\delta t}}\right)}^T\mathrm{.}$$

  
• Next the problem is formulated using this background knowledge.
3 Problem formulation
• A total of E objects are observed within the set

  

  . The vectors oⁱ = (xⁱ ,yⁱ ) ^T and θ ⁱ are the i^th object's position- and rotation-vector respectively. Tⁱ represents the i^th object's rotation-matrix dependent on θ ⁱ . The dynamical process of object i with state-vector sⁱ , process-noise wⁱ and measurement-vector mⁱ is defined with the following state-space model: $$s_k^i:=A_{{\tau }_k}{s}_{k-1}^i+B_{{\tau }_k}{w}_{k-1}^i,$$

  

  $$m_k^i:=\left(\begin{array}{c}{o}_k^i\\ {\theta }_k^i\end{array}\right)=C{s}_k^i+v_k^i,$$

  

  with $$C:=\left(I^{3\times 3}0\right),p\left(w_k^i\right):=G\left(w_k^{i}0Q^i\right),p\left(v_k\right):=f\left(v_k\right)\mathrm{.}$$

  
• The definition of the elements of state-vector sⁱ (t), also shown in Figure 18, are: $$s^i\left(t\right):={\left(x^i\left(t\right)y^i\left(t\right){\theta }^i\left(t\right)\frac{d{x}^i\left(t\right)}{\mathit{dt}}\frac{d{y}^i\left(t\right)}{\mathit{dt}}\right)}^T\mathrm{.}$$

  
• The objects are observed in

  

  by a camera or a network of sensors. For that M 'detection' points are marked within

  

  and collected in the set D ⊂

  

  . The position of a detection point is denoted as d ∈ D. The k^th detection of the system generates the observation vector $z_k^i\in \left\{{ℝ}_nℝ\right\}$

  

  if the edge of the i^th object covers one of the detection points d_k ∈ D at time t_k : $$z_k^i:=\left(d_k{t}_k\right),{}_{}\mathit{if}{\left(T_k^i\right)}^{-1}\left(d_k-o_k^i\right)\in S_E^i\mathrm{.}$$

  
• However, the system does not know which object was detected for it can be any object. As a result the system will not generate $z_k^i$

  

  but a general observation vector z_k ∈ {

  

  ,

  

  }, which is yet to be associated with an object. Therefore, due to the k^th detection, the observation vector z_k is generated whenever one of the E object covers a detection point d_k ∈ D at time t_k : $$z_k:=\left(d_k{t}_k\right),{}_{}\mathit{if}\left\{\exists i:{\left(T_k^i\right)}^{-1}\left(d_k-o_k^i\right)\in S_E^i\right\}\mathrm{.}$$

  
• From equations (9) and (10) we conclude that $z_k^i$

  

  of (9) is the result after the received observation vector z_k (10) is associated with object i. Notice that both definitions of z_k and $z_k^i$

  

  assume that the object's edge is detected exactly at a detection point d. In reality the detection will be affected by noise. The object therefore has some probability to be detected at a position υ ∈

  

  which is close to d. This is modeled by defining that the object's position at the instant of the detection, i.e. υ, is a random vector with mean d and covariance ε ∈

  

  : $$p\left(\upsilon \right):=G\left(\upsilon d\mathit{\epsilon I}\right)\mathrm{.}$$

  
• Figure 18 shows an example of object i which is detected by multiple detection points. The covariance ε of each detection point is also indicated.
• The sampling method of the observation vectors z _0:k is a form of event sampling [4, 5,7]. For a new observation vector is sampled whenever an event, i.e. object detection, takes place. With these event samples all N objects are to be tracked. To accomplish that three methods are needed. The first one is the association of the new observation-vector z_k to an object i and therefore denote it with $z_k^i\mathrm{.}$

  

  Suppose that all associated observation-vectors $z_n^i$

  

  are collected in the set $Z_k^i\in \left\{z_{0:k}\right\}\mathrm{.}$

  

  Then the second method is to estimate $m_k^i$

  

  from the observation-set $Z_k^i\mathrm{.}$

  

  This is used in the third method, which is a state-estimator.
• Present association methods are: Gating and Nearest Neighbor (NN) [2], (Joint) Probabilistic Data Association ((J)DPA) [2,8], Multiple Hypothesis Tracker (MHT) [9] and the MCMCDA [10]. Although these can be transformed for associating the event samples z _0:k , this paper will show that the estimation of $m_k^i$

  

  results in a probability that z_k is in fact $z_k^i,$

  

  i.e, $\mathit{Pr}\left(z_k=z_k^i\right)\mathrm{.}$

  

  Therefore the problem which is covered in this paper is the estimation of $m_k^i$

  

  from the set $Z_k^i,$

  

  which also results in the probability $\mathit{Pr}\left(z_k=z_k^i\right)\mathrm{.}$

  

  For that we assume that the shape of the object is known and that it is samples as shown in Section 2.3. The state-estimation is not covered in this paper, although it is used in the application example. Before going into the mathematical details of the estimation we will first describe the approach that is taken.
4 Approach measurement estimation
• In the problem formulation we stated that $Z_k^i$

  

  is defined as the set with all observation-vectors from z _0:k that were associated with object i. We will first redefine this set before continuing with the approach for estimating $m_k^i\mathrm{.}$

  
• The set $Z_k^i$

  

  is defined as the set of all observation-vectors z_n which were associated with object i, from which their detection point is still covered by the object. We will first show how this is done. At time step k we have the observation-set $Z_{k-1}^i$

  

  and the observation z_k was associated to object i, i.e. $z_k^i\mathrm{.}$

  

  Now if the object's edge is detected at d_k for the first time, then $z_k^i$

  

  is added to the set $Z_{k-1}^i\mathrm{.}$

  

  However, if the object's edge is detected at d_k for the second time, then $z_k^i\mathrm{.}$

  

  is not added to the set $Z_{k-1}^i$

  

  and the vector $z_n^i,$

  

  for which holds that d_n = d_k. is removed from $Z_{k-1}^i\mathrm{.}$

  

  This because in the second case, it means that object i drove off the detection point positioned at d_n = d_k . Therefore $Z_k^i$

  

  is defined as: $$Z_k^i:=\{\begin{array}{ll}{Z}_{k-1}^i\cap z_k,& \mathit{if}{d}_k\ne d_n,\forall z_n^i\in Z_{k-1}^i,\\ Z_{k-1}^i\backslash z_n,& \mathit{if}\left\{\exists n|d_k=d_n\wedge z_n^i=\left(d_n{t}_n\right)\in Z_{k-1}^i\right\}\mathrm{.}\end{array}$$

  
• With this definition of $Z_k^i$

  

  the approach for estimating $m_k^i,$

  

  i.e. $\left(m_k^i|Z_k^i\right),$

  

  is given. For clarity we assume that the object's shape is rectangular and that all its detection points are denoted with d_n , with n ∈ N ⊂ [0,k].
1. 1. The first step is to position the object on each detection point d_n and mirror its set S into the set O_n , as shown in Figure 19 for a single detection. This way we transform the points that are covered by the object, into possible vectors of the object's position $o_k^i\in O_n$
   
   given that it is detected at the detection point d_n .
2. 2. The second step, graphically depicted in Figure 20, is to turn all sets O_n simultaneously around their detection point d_n . This way, each possible orientation ${\theta }_k^i$
   
   of the object results in a corresponding possible object's position $o_k^i\mathrm{.}$
   
   For $o_k^i$
   
   must be inside all the sets O_n , ∀n ∈ N, and therefore thus inside the intersection of all sets O_n , ∀n ∈ N, which is denoted as O_N.
• Therefore if we apply these two steps for a number of orientations ${\theta }_k^i,$

  

  then at each orientation we have a set O_N which has to contain the object's position $o_k^i\mathrm{.}$

  

  From all these orientations we can calculate $p\left(m_k^i|Z_k^i\right)$

  

  as shown in the next section.
5 Measurement estimation
• Estimation of the measurement-vector $m_k^i$

  

  given the observation set $Z_k^i$

  

  results in calculating $p\left(m_k^i|Z_k^i\right)\mathrm{.}$

  

  Because both $m_k^i$

  

  and $Z_k^i$

  

  always belong to the same object and at sample instant k throughout this section we will remove the sub- and superscripts i and k in the rest of this section. Therefore we have; $m_k^i\to m$

  

  and $Z_k^i\to Z\mathrm{.}$

  

  The set Z consists of the observation vectors z_n , for all n ∈ N ⊂ [0,k], that were associated to the same object.
• Although the measurement vector is defined as m = (o, θ) ^T , with o ∈

  

  and θ ∈

  

  , the detection point at time-step n are defined as d_n ∈

  

  . Meaning that the objects orientation is not directly. However, because every observation vector z_n ∈ Z detects the object for one and the same θ, the PDF p(m|Z) is approximated by sampling in θ, i.e.: $$p\left(m|Z\right)\approx {\left({\displaystyle \sum _{l=1}^L}{\alpha }_l\right)}^{-1}{\displaystyle \sum _{l=1}^L}{\alpha }_{l}p\left(o|Z,\theta =\mathit{l\Delta }\right),$$

  

  $${\mathrm{with}}_{}\mathit{\Delta }:=\frac{2\pi }{L}\mathrm{and}{\alpha }_l:=\mathit{Pr}\left(\theta =\mathit{l\Delta }|Z\right)\mathrm{.}$$

  
• The main aspect of equation (13a) is to determine p(o|Z,θ). To do that we define the set O_n (θ) ∈

  

  to be equal to all possible object positions o, given that the object is detected at position d_n ∈ z_n (∈ Z) and that the object's rotation is equal to θ. The determination of O_n (θ) ∈

  

  is presented in the n the next section. Therefore, if one object is detected at multiple detection points d_n , ∀n ∈ N, then the set of all possible object positions o given a certain θ equals O_n (θ): $$O_N\left(\theta \right):=\underset{n\in N}{\cap }{O}_n\left(\theta \right)\mathrm{.}$$

  
• Equation (14) is graphically explained in Figure 21 for two different values of θ and N = {1,2}.
• Both p(o|Z,θ) and α _l are related to the set O_N (θ) due to the fact that it O_N (theta) defines the set of possible object positions o for a given θ. To calculate p(o|Z, θ) and α _l we define the functions f(o|Z,θ) and g(o|Z,θ): $$\begin{array}{ll}g\left(o|z_n,\theta \right)& :=\{\begin{array}{cc}0\mathit{if}& o\notin O_N\left(\theta \right),\\ 1\mathit{if}& o\in O_N\left(\theta \right),\end{array}\\ g\left(o|Z,\theta \right)& :={\mathrm{\Pi }}_{n\in N}f\left(o|z_n,\theta \right)=\{\begin{array}{cc}0\mathit{if}& o\notin O_N\left(\theta \right),\\ 1\mathit{if}& o\in O_N\left(\theta \right),\end{array}\end{array}$$

  
• Therefore the PDF p(o|Z,θ) and probability α _l are: $$p\left(o|Z,\theta \right):=\frac{g\left(o|Z,\theta \right)}{{\int }_{-\infty }^{\infty }g\left(o|Z,\theta \right)do},{\alpha }_l:=\frac{{\int }_{-\infty }^{\infty }g\left(o|Z,\theta \right)do}{{\int }_{-\infty }^{\infty }f\left(o|z_n,\theta \right)do}\mathrm{.}$$

  
• With (16) both p(m|Z) is calculated according to (13). The rest of this section is divided into two parts. The first part derives the probability function based on a single detection, i.e. f(o|z_n ,θ). While the second part derives the probability function based on a multiple detections, i.e. g(o|Z,θ).
5.1 Single event detection
• In order to derive f(o|z_n ,θ) we will use the set Λ, defined in 2.3, which contains the sampled positions λ _i that are covered by the object if o = 0 and θ = 0. Notice that if the object covers the origin, i.e. (x,y) ^T = 0, then the possible values of the object position o are given by the set -Λ. This is graphically depicted in Figure 22 (left). From that we can conclude that if the object covers the detection point d_n , given a certain orientation θ and rotation-matrix T, the sampled set Λ can be transformed into a sampled set of O_n , denoted with Ô _n : $${\hat{O}}_n:=\left[{\hat{o}}_1{\hat{o}}_2\dots {\hat{o}}_K\right],\mathrm{with}{\hat{o}}_i:=d-T{\lambda }_i\in O_n\mathrm{.}$$

  
• Figure 22 (right) graphically depicts the determination of Ô _n from the set Λ for a given θ and detection point d_n .
• The function f(o|z_n ,θ), as defined in (15), is approximated by placing a Gaussian function at each sampled position ô _i ∈ Ô _n with a certain covariance dependent on the grid-size r: $$f\left(o|z_n,\theta \right)\approx 2\pi {\gamma }^2{\displaystyle \sum _{i=1}^K}G\left(o{\hat{o}}_i{\mathit{\gamma I}}^{2\times 2}\right),$$

  

  with, $$\gamma :={\left(\frac{2\mathrm{<figure-callout\; id="2"\; label="\Delta \; K"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">\Delta </figure-callout>}}{K}\right)}^2\left(0.25-0.05e^{-\frac{4\left(K-1\right)}{15}}-0.08e^{-\frac{4\left(K-1\right)}{180}}\right)\mathrm{.}$$

  
• The approximation of (18) assumes that the object is detected exactly at d_n . In Section 3 we stated that the detection can be a bit of a detection point. The PDF that the object is detected at position υ ∈

  

  given the detection point d_n is defined in (11). Inserting this uncertainty into (18) results in the final f(o|z_n ,θ): $$f\left(o|z_n,\theta \right)\approx 2\pi {\gamma }^2{\displaystyle \sum _{i=1}^K}{\int }_{-\infty }^{\infty }G\left(\upsilon d\mathit{\epsilon I}\right)G\left(o,\upsilon -T{\lambda }_{i,\mathit{\gamma I}}\right)\mathit{d\upsilon }\mathrm{.}$$

  
• Equation (19) is solved with the following Proposition and the fact that G(x,a + b,C) = G(x-b,a,C):
• Proposition 1. Let there exist two Gaussian functions of the random vectors x ∈

  

  and m ∈

  

  and the matrix Γ ∈

  

  ; G(x,u,U) and G(m,Γx,M). Then they have the following property: $${\int }_{-\infty }^{\infty }G\left(xuU\right)G\left(m\mathit{\Gamma x}M\right)dx=G\left(\mathit{\Gamma u}\mathit{,}\mathit{m}\mathit{,}{\mathit{\Gamma U\Gamma }}^{\mathit{T}}+M\right)\mathrm{.}$$

  

  Proof. The proof can be found in Section 9.
• Applying Proposition 1 to (19) results in: $$f\left(o|z_n,\theta \right)\approx 2\pi {\gamma }^2{\displaystyle \sum _{i=1}^K}G\left(o{\hat{o}}_{i}R\right),{\mathrm{with}}_{}R:=\left(\epsilon +\gamma \right)I^{2\times 2}\mathrm{.}$$

  
• From f(o|z _n,θ) based on a single detection, the next step to multiple detections, i.e. g(olZ,θ), is taken.
5.2 Multiple event detections
• The aim of this section is to calculate the function g(o|Z,θ) by substituting equation (21) in the definition of g(o|Z, θ) as shown in (15): $$g\left(o|Z,\theta \right)\approx {\displaystyle \prod _{n\in \mathrm{<figure-callout\; id="2"\; label="N"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">N</figure-callout>}}}2\pi {\gamma }^2{\displaystyle \sum _{i=1}^K}G\left(o{\hat{o}}_i^{n}R\right),\mathrm{with}{\hat{o}}_i^n:=d_n-T{\lambda }_i\mathrm{.}$$

  
• If N contains m elements, then calculating equation (22) would result in K^m products of m Gaussian functions and sum them afterwards. This would take too much processing power if m is large. That is why equation (22) is calculated differently.
• Instead of using all detection points d_n we will use a subset of them. The derivation of this subset is graphically depicted in Figure 23 for N = {1,2}. For that consider the rectangular set

  

  ∈

  

  of Section 2.3 defined by its corners [c ₁, c ₂, c ₃, c ₄]. For each detection point d_n we define the set

  

  (θ) ⊂

  

  with corner-points [c _{n 1} (θ), c _{n 2} (θ), c _{n 3} (θ), c _{n 4} (θ)] defined as: $$c_{n_i}\left(\theta \right):={\mathit{Tc}}_i+d_n\mathrm{.}$$

  
• Let us define the rectangular set

  

  (θ) ⊂

  

  as the intersection of the sets

  

  (θ), ∀n ∈ N, i.e.: $${ℂ}_N\left(\theta \right):=\underset{n\in N}{\cap }{ℂ}_N\left(\theta \right)\mathrm{.}$$

  
• Meaning that each detection point d_n defines a rectangular set denoted with

  

  (θ) dependent on rotation θ. The intersection of all these rectangular sets is defined with the set

  

  (θ).
• In the beginning of this section we defined two different sets shown in Figures 19 and 20. The first set, O_n (θ), shown in Figure 19 defines all possible objet positions o based on a single detection at d_n. The second set, i.e. O_N (θ), shown in Figure 20, defines all possible object positions o based on all detections at d_n, ∀n ∈ N. Notice that as a result O_n (θ) ⊂

  

  (θ) and O_N (θ) ⊂

  

  (θ). Meaning that only within the set

  

  (θ) all the functions f(o|z_n, θ) have an overlapping area in which they are 1. Outside

  

  (θ) there is always at least one f(o|z_n, θ) which is 0 and therefore makes g(o|Z, θ) outside

  

  (θ) equal to 0. Therefore g(o|Z, θ) of (22) can be approximated by taking only those Gaussians of the functions f(o|z_n, θ) into account of which their mean, i.e. ${\hat{o}}_i^n,$

  

  is close or in the set

  

  (θ). We define that close to

  

  (θ) means a distance of at most γ + ε, which defined R in (21). The function g(o|Z, θ) of (22) is therefore approximated as: $$g\left(o|Z,\theta \right)\approx {\displaystyle \prod _{n\in \mathrm{<figure-callout\; id="2"\; label="N"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">N</figure-callout>}}}2\pi {\gamma }^2{\displaystyle \sum _{i\in N}^K}G{\left(o{\hat{o}}_i^{n}R\right)}_{†}$$

  

  with $$\left\{{\mathrm{I}}^n:\left|\left({\hat{o}}_i^n-{ℂ}_N\left(\theta \right)\right)\right|\le \left(\gamma +\epsilon \right),\forall i\in {\mathrm{I}}^j\right\}\mathrm{.}$$

  

  We can even decrease the number of Gaussians of (25) even further. This because if for a certain detection point d_n it holds

  

  (θ) ⊂ O_n (θ), it means that when we remove the detection point d_n it will not affect the set

  

  (θ). Therefore equation (25) is reduced to: $$g\left(o|Z,\theta \right)\approx {\displaystyle \prod _{n\in \mathrm{<figure-callout\; id="2"\; label="N"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">N</figure-callout>}}}2\pi {\gamma }^2{\displaystyle \sum _{i\in {\mathrm{I}}^n}^{}}G\left(o{\hat{o}}_i^{n}R\right),$$

  

  with $$\left\{\mathbb{N}\subset N:{ℂ}_N\left(\theta \right)\subset O_n\left(\theta \right),\forall n\in \left\{N\backslash \mathbb{N}\right\}\right\}\mathrm{.}$$

  
• The calculation of (26) is done by applying the following two propositions. The first one, i.e. Proposition 2, shows how to rewrite a product of a summation of Gaussians into a summation of a product of Gaussians. The second one, i.e. Proposition 3, proofs that a product of Gaussians results in a single Gaussian.
• Proposition 2. The product of a summation of Gaussians can be written into a summation of a product of Gaussian: $${\displaystyle \prod _{j=1}^C}{\rho }_j{\displaystyle \sum _{i=1}^{C_j}}G\left(x{x}_j^{i}R\right)=\left({\displaystyle \prod _{j=1}^C}{\rho }_j{C}_j\right){\displaystyle \sum _{m=1}^{\prod _{j=1}^C{C}_j}}{\displaystyle \prod _{j=1}^C}G\left(x{x}_j^{f\left(jm\right)}R\right),$$

  

  $$\begin{array}{l}f\left(jm\right):=\lceil \frac{m-\sum _{r=1}^j\left(f\left(r-1,m\right)-1\right)\left(\prod _r^{C+1}{\rho }_r\right)}{\prod _{r=j+1}^{C+1}{\rho }_r}\rceil ,\\ \mathrm{<figure-callout\; id="0"\; label="f"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">f</figure-callout>}\left(0m\right):=1,{\rho }_{j+1}:=1.\end{array}$$

  
• The proof is given by writing out the left hand side of (27a) and restructuring it.
• Proposition 3. The product of Gaussians is again a Gaussian: $${\displaystyle \prod _{j=1}^C}G\left(x{x}_{j}R\right)=\mathit{\beta G}\left(\frac{x,\sum _{j=1}^C{x}_j}{C}\frac{R}{C}\right)\mathit{and\; \beta }={\displaystyle \prod _{j=2}^C}G\left(x_j\frac{\sum _{n=1}^{j-1}{x}_n}{j-1}\frac{\mathit{jR}}{j-1}\right)\mathrm{.}$$

  
• The proof is given in Section 10.
• Now applying Propositions 2 and 3 on (26) results in a solution of g(o|Z, θ) as a summation of Gaussians of the form: $$g\left(o|Z,\theta \right)={\displaystyle \sum _{i=1}^H}{\beta }^i\left(\theta \right)G\left(o,o^i\left(\theta \right),R^i\left(\theta \right)\right),$$

  
• Equation (29) is approximated as a single Gaussian function: $$g\left(o|Z,\theta \right)\approx \stackrel{}{\beta }\left(\theta \right)G\left(o,\stackrel{}{o}\left(\theta \right),\stackrel{}{R}\left(\theta \right)\right),$$

  

  $$\begin{array}{l}\stackrel{}{\beta }\left(\theta \right):={\displaystyle \sum _{i=1}^H}{\beta }^i\left(\theta \right),\stackrel{}{o}\left(\theta \right):={\displaystyle \sum _{i=1}^H}\frac{{\beta }^i\left(\theta \right)}{\stackrel{}{\beta }}{o}^i\left(\theta \right),\\ \stackrel{}{R}\left(\theta \right):={\displaystyle \sum _{i=1}^H}\frac{{\beta }^i\left(\theta \right)}{\stackrel{}{\beta }}\left(R^i\left(\theta \right)+\left(\stackrel{}{o}\left(\theta \right)-o^i\left(\theta \right)\right){\left(\stackrel{}{o}\left(\theta \right)-o^i\left(\theta \right)\right)}^T\right)\mathrm{.}\end{array}$$

  
• With the result of (30) we can approximate g(o|Z, θ). In order to calculate the PDF p(m|Z), equation (30) is substituted into equation (16) together with f(o|Z_n, θ) of (21) to calculate p(o|Z, θ) and α₁. Substituted these results into (13) gives: $$p\left(m|Z\right)={\displaystyle \sum _{l=1}^L}\frac{\stackrel{}{\beta }\left(\mathit{l\Delta }\right)}{\sum _{l=1}^L\stackrel{}{\beta }\left(\mathit{l\Delta }\right)}G\left(o,\stackrel{}{o}\left(\mathit{l\Delta }\right),\stackrel{}{R}\left(\mathit{l\Delta }\right)\right)\mathrm{.}$$

  
• As was mentioned in the problem formulation, the PDF p(m|Z) also gives us the probability that a new observation vector is generated by an certain object i. This is discussed in the next section.
6 Detection association
• The total probability that a new observation vector z_k is generated by object i is equal to the total probability of the measurement-vector $m_k^i$

  

  given the observation set $\left[Z_{k-1}^i{z}_k\right]\mathrm{.}$

  

  For this probability we can use $p\left(m_k^i|Z_{k-1}^i,z_k\right)$

  

  which is equal to equation (31). The definition of a PDF is that its total probability, i.e. its integral from -∞ to ∞, is equal to 1. To make sure that $p\left(m_k^i|Z_{k-1}^i,z_k\right)$

  

  of equation (31) has a total probability of 1, it is divided by its true probability $\mathit{Pr}\left(m_k^i|Z_{k-1}^i,z_k\right)\mathrm{.}$

  

  In order to be able to compare these different measurement-vector per object, we normalize each probability with the surface covered by the object. As a result, $\mathit{Pr}\left(z_k=z_k^i\right)$

  

  is equal to: $$\mathit{Pr}\left(z_k=z_k^i\right)=\frac{1}{2\pi {\left({\gamma }^i\right)}^2{K}^i}{\displaystyle \sum _{l=1}^L}\stackrel{}{\beta }\left(\mathit{l\Delta }\right)\mathrm{.}$$

  
• The variables γⁱ and Kⁱ are equal to γ and K respectively, which define the approximation of the function $f(m_k^i|z_n^i,{\theta }_k^i$

  

  as shown in (18). With the probability of (32) one can design a method which associates an observation-vector due to a new detection, to its most probable object i. Although the estimation method requires a certain amount of processing power, one can reduce this by reducing the number of samples in the set A. Meaning that association and estimation can be done with different sizes of A. Moreover, if the objects have a rectangular shape, then with some tricks one can reduce the amount of processing power to a level at which both association as well as estimation can run real-time.
• Now that both the measurement estimation as well as the detection association are designed, both are tested in a multiple object tracking application.
7 Application example
• As an application example we take a parking lot of 50 by 50 meters with a network of wireless sensors distributed randomly along the road's surface. Each sensor can detect a crossing vehicle. A total of 2500 sensors was used resulting in a density of one sensor per square meter. The vehicles are all assumed rectangular shaped objects with a length of 5 meters and a width of 2 meters. A total of 4 vehicles manoeuvre within the parking lot and are tracked using a data-associator followed with an event state-estimator.
• The simulation case is made such that it contains two interesting situation. One in which two vehicles cross each other in parallel and one where two vehicles cross perpendicular. For comparison the objects are tracked using two different association methods. The first one is a combination of Gating and detection association of 6. The second one is a combination of Gating and Nearest Neighbor.
• The result of the detection associator (DA) for both crossings is shown in Figure 24 while the result of the Nearest Neighbor (NN) associator is shown in Figure 25. In both results the real object is plotted in a thick, solid line while its estimated one is plotted in a thin, solid line. The associated detections of each object are given with a symbol which is different for each object; '□' if associated with vehicle 1, '□' if associated with vehicle 2, '∇' if associated with vehicle 3 and '*' if associated with vehicle 4. Figure 24 shows with the DA all detections were correctly associated to the one object, while If NN is used as an association method, we see that a lot of incorrect associated detections. Therefore we can concluded that using the detection association of 6 results in less estimation-error compared to NN.
• A second simulation is done to compare the percentage of incorrect associated detections. Again for the both DA as well as NN only now 4 different amount of detection points were used: 3000, 2500, 2000 and 1500. This table shows that the detection association has a better performance compared to Nearest Neighbor. Table 1. Percentage of incorrect association

  amount of detection points	DA	NN
  3000	0%	4.5%
  2500	0%	5.6%
  2000	0%	7.8%
  1500	0%	2.2%

8 Conclusions
• This paper presents a method for estimating the position- and rotation-vectar of objects from spatially, distributed detections of that object. Each detection is generated at the event that the edge of an object crosses a detection point. From the estimation method a detection associator is also designed. This association method calculates the probability that a new detection was generated by an object i.
• An example of a parking lot shows that the detection association method has no incorrect associated detections in the case that two vehicles cross each other both in parallel as well as orthoganal. If the association method of Nearest Neighbor was used, a large amount of incorrect associated detections were noticed, resulting in a higher state-estimation error.
• The data-assimilation can be further improved with two adjustments. The first one is replacing the set S with S_E only at the time-instants that the observation vector is received. The second improvement is to take the detection points that have not detected anything also in account.
References
• 1. A. Poore and S. Gadaleta, "Some assignment problems arising from multiple target tracking," Mathematical and Computer Modelling, vol. 43, pp. 1074-1091, 2006.
• 2. Y Bar-Shalom and R. Li, Multitarget-Multisensor Tracking: Principles and Techniques. YBS, 1995.
• 3. S. Blackman and R. Popoli, Design and Analysis of Modern Tracking Systems. Artech House, Norwood, MA, 1999.
• 4. V. Nguyen and Y Suh, "Improving estimation performance in Networked Control Systems applying the Send-on-delta transmission method," Sensors, vol. 7, pp. 2128-2138, 2007.
• 5. K. Åström and B. Bemhardsson, "Comparison of Riemann and Lebesque sampling for first order stochastic systems," in 41st IEEE Conference on Decision and Control, Las Vegas, Nevada, USA, 2002.
• 6. N. Johnson, S. Kotz, and A. Kemp, Univariate discrete distributions. JOHN WILEY and SONS, Inc, 1992.
• 7. M. Miskowicz, "Send-on-deltn concept: an event-based data-reporting strategy," Sensors, vol. 6, pp. 49-63, 2006.
• 8. R. Karlsson and F. Gustafsson, "Monte Carlo data association for multiple target tracking," in IEE International Seminar on Target Tracking: Algorithms and Applications, 2001.
• 9. B. Ristic, S. mlampalam, and N. Gordon, Beyond the Kahnan filter: Particle filter for tracking applications, 2002.
• 10. O. Songhwai, S. Sastry, and L. Schenato, "A Hierarchical Multiple-Target Tracking Algorithm for Sensor Networks," in Proc. of the 2005 IEEE International Conference on Robotics and Automation, Barcelona, Spain, 2005.
• 11. K. Mardia, J. Kent, and 1. Bibby, Mutlivariate analysis. Academic press, Inc. London, 1979.
• 12. D. Kelly, Introduction to probability. Macmillan Publishing Company, 1994.
• 13. H. Durant-Whyte, B. Rao, and H. Hu, "Towards a fully decentralized architecture for multisensor data fusion," in 1990 IEEE Int. Conf. on Robotics and Automation, Cincinnati, Ohio, USA, 1990, pp. 1331-1336.
9 Proof of Proposition 1
• Proof. Defined are two Gaussian functions with the vectors x ∈

  

  , u ∈

  

  , m ∈

  

  and matrices U ∈

  

  , M ∈

  

  , Γ ∈

  

  : G(x,u,U) and G(m,Γx,M). Suppose we define the following PDFs and relation of m with some c ∈

  

  : $$m=\mathit{\Gamma x}+c,$$

  

  $$p\left(c\right):=G\left(m0M\right)\mathrm{and}p\left(x\right):=G\left(xuU\right)\mathrm{.}$$

  
• Then from probability theory [6] p(m) is equal to: $$p\left(m\right):={\int }_{-\infty }^{\infty }p\left(m|x\right)p\left(x\right)dx$$

  

  $$={\int }_{-\infty }^{\infty }G\left(m\mathit{\Gamma x}M\right)G\left(xuU\right)dx\mathrm{.}$$

  
• Applying theorem 3.2.1 of [11] on (33b) we have that p(Γx) = G(Γx,Γu,ΓUΓ ^T ). Now if we have the random vectors a ∈

  

  and b ∈

  

  with p(a) = G(a ₁ u ₁ ,U ₁) and p(b) = G(b,u ₂,U ₂) then they have the property p(a+b) = G(a+b,u ₁+u ₂,U ₁+U ₂) as proven in [12]. Applying this on (33a) results in: $$p\left(m\right)=G\left(m\mathit{,}\mathit{\Gamma u}\mathit{,}{\mathit{\Gamma U\Gamma }}^{\mathit{T}}+M\right),$$

  

  $$\Rightarrow {\int }_{-\infty }^{\infty }G\left(m\mathrm{\Gamma x}M\right)G\left(xuU\right)dx=G\left(m\mathit{,}\mathit{\Gamma u}\mathit{,}{\mathit{\Gamma U\Gamma }}^{\mathit{T}}+M\right)\mathrm{.}$$

  
10 Proof of Proposition 3
• Proof. A product of Gaussians can be written as: $${\displaystyle \prod _{j=1}^N}G\left(x{x}_{j}R\right)=G\left(x{x}_{N}R\right){\displaystyle \prod _{j=1}^{N-1}}G\left(x{x}_{j}R\right),$$

  

  $$={\beta }_{N-1}G\left(x{x}_{N}R\right)G\left(x\frac{\sum _{j=1}^{N-1}{x}_j}{N-1}\frac{R}{N-1}\right)\mathrm{.}$$

  
• From Proposition 1 and the Kalman filter in Information form [13], a product of 2 Gaussians equals: $$G\left(xuU\right)G\left(mxM\right)=G\left(m,u,U+M\right)G\left(xdD\right),$$

  

  $$\mathrm{with}{}_{}{D}^{-1}=U^{-1}+M^{-1},{}_{}d={\mathit{DU}}^{-1}u+{\mathit{DM}}^{-1}m\mathrm{.}$$

  
• Applying (37) on (36b), together with the fact that G(x,y,Z) = G(y,x,Z) we have: $${\displaystyle \prod _{j=1}^N}G\left(x{x}_{j}R\right)={\beta }_{N-1}G\left(x_N\frac{\sum _{j=1}^{N-1}{x}_j}{N-1}\frac{\mathit{NR}}{N-1}\right)G\left(x\frac{\sum _{j=1}^N{x}_j}{N}\frac{R}{N}\right)\mathrm{.}$$

  
• Equation (38) is equal to (28) for: $${\beta }_N={\beta }_{N-1}G\left(x_N\frac{\sum _{j=1}^{N-1}{x}_j}{N-1}\frac{\mathit{NR}}{N-1}\right),$$

  

  $$={\displaystyle \prod _{i=2}^N}G\left(x_i\frac{\sum _{j=1}^{N-1}{x}_j}{i-1}\frac{\mathit{iR}}{i-1}\right)\mathrm{.}$$

  
A2: On Event Based State Estimation
• Summary. To reduce the amount of data transfer in networked control systems and wireless sensor networks, measurements are usually taken only when an event occurs, rather that at each synchronous sampling instant. However, this complicates estimation and control problems considerably. The goal of this paper is to develop a state estimation algorithm that can successfully cope with event based measurements. Firstly, we propose a general methodology for defining event based sampling. Secondly, we develop a state estimator with a hybrid update, i.e. when an event occurs the estimated state is updated using measurements; otherwise the update is based on the knowledge that the monitored variable is within a bounded set used to define the event. A sum of Gaussians approach is employed to obtain a computationally tractable algorithm.
1 Introduction
• Different methods for state estimation have been introduced during the last decades. Each method is specialized in the type of process, the type of noise or the type of system architecture. In this paper we focus on the design of a state estimator that can efficiently cope with event based sampling. By even sampling we mean that measurements are generated only when an a priori defined event occurs in the data monitored by sensors. Such an effective estimator is very much needed in both networked control systems and wireless sensor networks (WSNs) [1]. Especially in WSNs, where the limiting resource is energy, data transfer and processing power must be minimized. The existing estimators that could be used in this framework are discussed in Section 4. For related research on event based control, the interested reader is referred to the recent works [2], [3].
• The contribution of this paper is twofold. Firstly, we introduce a general mathematical description of event based sampling. We assume that the estimator does not know when new measurements are available, which usually results in unbounded eigenvalues of its error-covariance matrix. To obtain an estimator with a bounded error-covariance matrix, we develop an estimation algorithm with hybrid update, which is the second main contribution. The developed event based estimator is updated both when an event occurs, with received measurements, as well as at sampling instants synchronous in time. Then the update is based on the knowledge that the monitored variable is within a bounded set used to define the event. In order to meet the requirements of a low processing power, the proposed state estimator is based on the Gaussian sum filter [4,5], which is known for its computational tractability.
2 Background notions and notation
• R defines the set of real numbers whereas the set

  

  defines the non-negative real numbers. The set

  

  defines the integer numbers and

  

  defines the set of non-negative integer numbers. The notation 0 is used to denote either the null-vector or the null-matrix. Its size will become clear from the context.
• A vector x(t) ∈

  

  is defined to depend on time t ∈

  

  and is sampled using some sampling method. Two different sampling methods are discussed. The first one is time sampling in which samples are generated whenever time t equals some predefined value. This is either synchronous in time or asynchronous. In the synchronous case the time between two samples is constant and defined as t_s ∈

  

  . If the time t at sampling instant k_a ∈

  

  is defined as t_ka , with t_0a := 0, we define: $$x_{k_a:=x\left(t_{k_a}\right)}\mathrm{and}{x}_{0_a:k_a}:=\left(x\left(t_{0_a}\right),x\left(t_{1_a}\right),\cdots ,x\left(t_{k_a}\right)\right)\mathrm{.}$$

  
• The second sampling method is event sampling, in which samples are taken when an event occurred. If t at event instant k_e ∈

  

  is defined as t_ke , with t _{0 e} , := 0, we define: $$x_{k_e:=x\left(t_{k_e}\right)}\mathrm{and}{x}_{0_e:k_e}:=\left(x\left(t_{0_e}\right),x\left(t_{1_e}\right),\cdots ,x\left(t_{k_e}\right)\right)\mathrm{.}$$

  
• A transition-matrix A _{t 2-t 1} ∈

  

  is defined to relate the vector u(t ₁) ∈

  

  to a vector x(t ₂) ∈

  

  as follows: x(t ₂) = A _{t 2 -t 1} u(t ₁).
• The transpose, inverse and determinant of a matrix A ∈

  

  are denoted as A^T, A ^-1 and |A| respectively. The i^th and maximum eigenvalue of a square matrix A are denoted as λ _i (A) and λ _max (A) respectively. Given that A ∈

  

  and B ∈

  

  are positive definite, denoted with A > 0 and B > 0, then A > B denotes A - B > 0. A ≥ 0 denotes A is positive semi-definite.
• The probability density function (PDF), as defined in [6] section B2, of the vector x ∈

  

  is denoted with p(x) and the conditional PDF of x given u ∈

  

  is denoted as p(x|u). The expectation and covariance of x are denoted as E[x] and cov(x) respectively. The conditional expectation of x given u is denoted as E[x|u]. The definitions of E[x], E[x|u] and cov(x) can be found in [6] sections B4 and B7.
• The Gaussian function (shortly noted as Gaussian) of vectors x ∈

  

  and u ∈

  

  and matrix P ∈

  

  is defined as G(x,u,P):

  

  ×

  

  ×

  

  →

  

  , i.e.: $$G\left(xuP\right)=\frac{1}{\sqrt{{\left(2\pi \right)}^n\left|P\right|}}{e}^{-0.5{\left(x-u\right)}^T{p}^{-1}\left(x-u\right)}\mathrm{.}$$

  
• If p(x) = G(x,u,P), then by definition it holds that E[x] = u and cov(x) = P.
• The element-wise Dirac-function of vector x ∈

  

  , denoted as δ(x):

  

  → {0,1}, satisfies: $$\delta \left(x\right)=\{\begin{array}{cc}0& \mathrm{if}x\not\equiv \underline{0},\\ 1& \mathrm{if}x\equiv \underline{0},\end{array}\mathrm{and}{\int }_{-\infty }^{\infty }\delta \left(x\right)dx=1.$$

  
• For a vectors x ∈

  

  and a bounded Borel set [7] Y ⊂

  

  , the set PDF is defined as Λ _Y (x):

  

  → {0, v} with v ∈

  

  defined as the Lebesque measure [8] of the set Y, i.e.: $${\mathit{\Lambda }}_{\mathit{\gamma }}\left(x\right)=\{\begin{array}{cc}0& {\mathrm{if}}_{}x\notin Y,\\ {\nu }^{-1}& {\mathrm{if}}_{}x\in Y\mathrm{.}\end{array}$$

  
3 Event sampling
• Many different methods for sampling a vector y(t) ∈

  

  can be found in literature. The one mostly used is time sampling in which the $k_n^m$

  

  sampling instant is defined at time t_ka := t _{ka -1} + τ _{ka -1} for some τ _{ka -1} ∈

  

  . Recall that if y(t) is sampled at t_a it is denoted as y_ka . This method is formalized by defining the observation vector z _{ka -1} := ${\left(y_{k_a-1}^T{t}_{k_a-1}\right)}^T\in {ℝ}^{q+1}$

  

  at sampling instant k _a-1. Let us define the set H_ka (z _{ka -1}) ⊂

  

  containing all the values that t can take between t _{ka -1} and t _{ka -1} + τ _{ka -1} i.e.: $$H_{k_a}\left(z_{k_a-1}\right):=\left\{t\in ℝ|t_{k_a-1}\le t<t_{k_a-1}+{\tau }_{k_a-1}\right\}\mathrm{.}$$

  
• Then time sampling defines that the next sampling instant, i.e. k_a , takes place whenever present time t exceeds the set H_ka (z _{ka -1}). Therefore z_ka is defined as: $$z_{k_a}:=\left(y_{k_a}{t}_{k_a}\right)\mathrm{if}t\notin H_{k_a}\left(z_{k_a-1}\right)\mathrm{.}$$

  
• In the case of synchronous time sampling τ _ka = t_s, ∀k_a ∈

  

  , which is graphically depicted in Figure 26. Notice that with time sampling, the present time t specifies when samples of y(t) are taken, but time t itself is independent of y(t). As a result y(t) in between the two samples can have any value within

  

  . Recently, asynchronous sampling methods have emerged, such as, for example "Send-on-Delta" [9,10] and "Integral sampling" [11]. Opposed to time sampling, these sampling methods are not controlled by time t, but by y(t) itself.
• Next, we present a general definition of event based sampling, which recovers the above mentioned asynchronous methods, for a particular choice of ingredients. Let us define the observation vector at sampling instant k_e - 1 as $z_{k_e-1}:={\left(y_{k_e-1}^T{t}_{k_e-1}\right)}^T\in$

  

  ${\mathrm{<figure-callout\; id="1"\; label="R\; q\; +"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">ℝ</figure-callout>}}^{q+}\mathrm{.}$

  

  With that we define the following bounded Borel set in time-measurement-space, i.e. H_ke (Z _{ke -1},t) ⊂

  

  , which depends on both z_ke -1 and t. In line with time sampling the next event instant, i.e. k_e , takes place whenever y(t) leaves the set H_ke (z _{ke -1}, t) as shown in Figure 27 for q = 2. Therefore z_ke is defined as: $$z_{k_e}:=\left(y_{k_e}{t}_{k_e}\right)\mathrm{if}y\left(t\right)\notin H_{k_e}\left(z_{k_e-1}t\right)\mathrm{.}$$

  

  The exact description of the set H_ke (z _{ke -1} , t) depends on the actual sampling method. As an example H_ke (z _{ke -1},t) is derived for the method "Send-on-Delta", with y(t) ∈

  

  . In this case the event instant k_e occurs whenever |y(t) - y _{ke -1}| exceeds a predefined level Δ, see Figure 28, which results in H_ke (z _{ke -1},t) = {y ∈

  

  | - Δ < y - y _{ke -1} < Δ}.
• In event sampling, a well designed H_ke (z _{ke -1}, t) should contain the set of all possible values that y(t) can take in between the event instants k_e - 1 and k_e . Meaning that if t _{ke -1} ≤ t < t_ke , then y(t) ∈ H_ke (z _{ke -1},t). A sufficient condition is that y _{k e-1} ∈ H_ke (z _{ke -1} ,t), which for "Send-an-Delta" results in y(t) ∈ [y _{ke -t} -Δ,y _e-1 +Δ] for all t_ke -1 ≤ t < t_ke .
4 Problem formulation: State estimation based on event sampling
• Assume a perturbed, dynamical system with state-vector x(t) ∈

  

  , process-noise w(t) ∈

  

  , measurement-vector y(t) ∈

  

  and measurement-noise v(t) ∈

  

  . This process is described by a state-space model with A _τ ∈

  

  , B _τ ∈

  

  and C ∈

  

  . An event sampling method is used to sample y(t). The model of this process becomes: $$x\left(t+\tau \right)=A_{\tau }x\left(t\right)+B_{\tau }w\left(t\right),$$

  

  $$y\left(t\right)=\mathit{Cx}\left(t\right)+v\left(t\right),$$

  

  $$z_{k_e}={\left(y_{k_e}^T{t}_{k_e}\right)}^T\mathrm{if}y\left(t\right)\notin H_{k_e}\left(z_{k_e-1}t\right),$$

  

  with $$p\left(w\left(t\right)\right):=G\left(w\left(t\right),0,Q\right)\mathrm{and}p\left(v\left(t\right)\right):=G\left(v\left(t\right),0,V\right)\mathrm{.}$$

  

  The state vector x(t) of this system is to be estimated from the observation vectors z _{0 e :ke} . Notice that the estimated states are usually required at all synchronous time samples k_a , with t_s = t _ka - t _{ka -1}, e.g., as input to a controller that runs synchronously in time. As such our goal is to construct an event-based state-estimator (EBSE) that provides an estimate of x(t) not only at the event instants t _ke but also at the sampling instants t_ka . Therefore, we define a new set of sampling instants t_n as the combination of sampling instants due to event sampling, i.e. k_e , and time sampling, i.e. k_a : $$\left\{t_{0:n-1}\right\}:=\left(t_{0_a:k_a-1}\right)\cup \left(t_{0_e:k_e-1}\right)\mathrm{and}{t}_n:=\{\begin{array}{cc}{t}_{k_a}& \mathrm{if}{t}_{k_a}<t_{k_e},\\ t_{k_e}& \mathrm{if}{t}_{k_a}\ge t_{k_e}\mathrm{.}\end{array}$$

  

  and $${\mathrm{<figure-callout\; id="0"\; label="t"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">t</figure-callout>}}_0<{\mathrm{<figure-callout\; id="1"\; label="t"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">t</figure-callout>}}_1<\cdots <t_n,x_n:=x\left(t_n\right),y_n:=y\left(t_n\right)\mathrm{.}$$

  
• The estimator calculates the PDF of the state-vector x_n given all the observations until t_n . This results in a hybrid state-estimator, for at time t_n an event can either occur or not, which further implies that measurement data is received or not, respectively. In both cases the estimated state must be updated (not predicted) with all information until t_n . Therefore, depending on t_n a different PDF must be calculated, i.e.: $$\mathrm{if}{t}_n=t_{k_a}\Rightarrow p\left(x_n|z_{0_e:k_e-1}\right)\mathrm{with}{t}_{k_e-1}<t_{k_a}<t_{k_e},$$

  

  $$\mathrm{if}{t}_n=t_{k_e}\Rightarrow p\left(x_n|z_{0_e:k_e}\right)\mathrm{.}$$

  

  The important parameters for the performance of any state-estimator are the expectation and error-covariance matrix of its calculated PDF Therefore, from (9) we define: $$x_{n|n}:=\{\begin{array}{cc}E\left[x_n|z_{0_e:k_e-1}\right]& \mathrm{if}{t}_n=t_{k_a}\\ E\left[x_n|z_{0_e:k_e}\right]& \mathrm{if}{t}_n=t_{k_e}\end{array}\mathrm{and}{P}_{n|n}:=\mathit{cov}\left(x_n-x_{n|n}\right)\mathrm{.}$$

  

  The PDFs of (9) can be described as the Gaussian G(x_n ,x _n|n , P_n |_n). The square root of the eigenvalues of P _n|n , i.e. $\sqrt{{\lambda }_i\left(P_{n|n}\right)},$

  

  define the shape of this Gaussian function. Together with x _n|n they indicate the bound which surrounds 63% of the possible values for x_n . This is graphically depicted in Figure 29 for the ID case and Figure 30 for a 2D case, in a top view. The smaller the eigenvalues λ _i (P _n|n ) are, the smaller the estimation-error is.
• As such, the problem of interest in this paper is to construct a state-estimator suitable for the general event sampling method introduced in Section 3 and which is computationally tractable. Furthermore, it is desirable to guarantee that P _n|n has bounded eigenvalues for all n.
• Existing state estimators can be divided into two categories. The first one contains estimators based on time sampling: the (a)synchronous Kalman filter [12, 13] (linear process, Gaussian PDF), the Particle filter [14] and the Gaussian sum filter [4, 5] (nonlinear process, non-Gaussian PDF). These estimators cannot be directly employed in event based sampling as if no new observation vector z_ke is received, then t_n - t_ke → ∞ and λ _i (P _{n|ke -1}) → ∞. The second category contains estimators based on event sampling. In fact, to the best of our knowledge, only the method proposed in [15] fits this category. However, this EBSE is only applicable in the case of "Send-on-Delta" event sampling and it requires that any PDF is approximated as a single Gaussian function. Moreover, the asymptotic property of P _n|n is not investigated in [15].
• In the next section we propose a novel event-based state-estimator, suitable for any event sampling method, together with a preliminary result on asymptotic analysis.
5 An event-based state estimator
• The EBSE estimates x_n given the received observation vectors until time t_n . Notice that due to the definition of event sampling we can extract information of all the measurement vectors y _{0: n}. For with t_i ∈ {t_0:n } and t_je ∈ {t_{0 e :ke} } it follows that: $$\{\begin{array}{cc}{y}_i\in H_{j_e}\left(z_{j_e-1}{t}_i\right)& \mathrm{if}{t}_{j_e-1}\le t_i<t_{j_e},\hfill \\ y_i=y_{j_e}& \mathrm{if}{t}_i=t_{j_e}\mathrm{.}\hfill \end{array}$$

  
• Therefore, from the observation vectors z _0e:ke and (11) the PDFs of the hybrid state-estimation of (9), with the bounded, Borel set Y_i ⊂

  

  , results in: $$p\left(x_n|{\mathrm{<figure-callout\; id="0"\; label="y"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">y</figure-callout>}}_0\in {\mathrm{<figure-callout\; id="0"\; label="Y"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">Y</figure-callout>}}_0,{\mathrm{<figure-callout\; id="1"\; label="y"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">y</figure-callout>}}_1\in {\mathrm{<figure-callout\; id="1"\; label="Y"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">Y</figure-callout>}}_1,\dots ,y_n\in Y_n\right)$$

  

  $$Y_i:=\{\begin{array}{cc}{H}_{j_e}\left(z_{j_e-1}{t}_i\right)& \mathrm{if}{t}_{j_e-1}<t_i<t_{j_e},\hfill \\ \left\{y_{j_e}\right\}& \mathrm{if}{t}_i=t_{j_e}\mathrm{.}\hfill \end{array}$$

  
• For brevity (12a) is denoted as p(x_n |y _0:n ∈ Y _0:n ) and with Bayes-rule [16] yields: $$p\left(x_n|y_{0:n}\in Y_{0:n}\right):=\frac{p\left(x_n|y_{0:n-1}\in Y_{0:n-1}\right)p\left(y_n\in Y_n|x_n\right)}{p\left(y_n\in Y_n|y_{0:n-1}\in Y_{0:n-1}\right)}\mathrm{.}$$

  
• To have an EBSE, with low processing demand, multivariate probability theory [17] is used to make (13) recursive: $$p\left(a|b\right):={\int }_{-\infty }^{\infty }p\left(a|c\right)p\left(c|b\right)dc\Rightarrow$$

  

  $$p\left(x_n|y_{0:n-1}\in Y_{0:n-1}\right)={\int }_{-\infty }^{\infty }p\left(x_n|x_{n-1}\right)p\left(x_{n-1}|y_{0:n-1}\in Y_{0:n-1}\right)d{x}_{n-1},$$

  

  $$p\left(y_n\in Y_n|y_{0:n-1}\in Y_{0:n-1}\right)={\int }_{-\infty }^{\infty }p\left(x_n|y_{0:n-1}\in Y_{0:n-1}\right)p\left(y_n\in Y_n|x_n\right)d{x}_n\mathrm{.}$$

  
• The calculation of p(x_n |y _0:n ∈ Y _0:n ) is done in three steps: 1. Assimilate p(y_n ∈ Y_n |x_n ) for both t_n = t_ke and t_n = t_ka . 2. Calculate p(x_n |y _0:n ∈ Y _0:n ) as a summation of N Gaussians. 3. Approximate p(x_n |y _0:n ∈ y_0:n ) as a single Gaussian function. The reason for this last step is to design an algorithm in which p(x _n|y_0:n ∈ Y_0:n ) is described by a finite set of Gaussians and therefore attain computational tractability. Notice that (13) gives a unified description of the hybrid state-estimator, which makes an asymptotic analysis of the EBSE possible, as it will be shown later in this section.
5.1 Step 1: measurement assimilation
• This section gives a unified formula of the PDF p(y_n ∈ Y_n |x_n ) valid for both t_n = t_ke and t_n = t_ka . From multivariate probability theory [17] and (7b) we have: $$p\left(y_n\in Y_n|x_n\right):={\int }_{-\infty }^{\infty }p\left(y_n|x_n\right)p\left(y_n\in Y_n\right)d{y}_n\mathrm{and}p\left(y_n|x_n\right)=G\left(y_n{\mathit{Cx}}_{n}V\right)\mathrm{.}$$

  
• The PDF p(y_n ∈ Y_n ) is modeled as a uniform distribution for all y_n ∈ Y_n . Therefore, depending on the type of instant, i.e. event or not, we have: $$p\left(y_n\in Y_n\right):=\{\begin{array}{cc}{\mathit{\Lambda }}_{H_{k_e}}\left(y_n\right)\hfill & \mathrm{if}{t}_{k_e-1}<t_n<t_{k_e},\hfill \\ \delta \left(y_n-y_{k_e}\right)\hfill & \mathrm{if}{t}_n=t_{k_e}\mathrm{.}\hfill \end{array}$$

  
• Substitution of (16) into (15) gives that p(y_n ∈ Y_n |x_n ) = G(y_ke ,Cx_n, V) if t_n = t_ke . However, if t_n = t_ka then p(y_n ∈ Y_n |x_n ) equals ${\mathrm{\Lambda }}_{{\mathit{H}}_{{\mathit{k}}_{\mathit{e}}}}\left({\mathit{y}}_{\mathit{n}}\right),$

  

  which is not necessarily Gaussian. Moreover, it depends on the set H_ke and therefore on the actual event sampling method that is employed. In order to have a unified expression of p(y_n ∈ Y_n |x_n ) for both types of t_n , independent of the event sampling method, ${\mathrm{\Lambda }}_{{\mathit{H}}_{{\mathit{k}}_{\mathit{e}}}}\left({\mathit{y}}_{\mathit{n}}\right)$

  

  can be approximated as a summation of N Gaussians, i.e. $$\begin{array}{c}{\mathit{\Lambda }}_{H_{k_e}}\left(y_n\right)\hfill \end{array}\approx {\displaystyle \sum _{i=1}^N}{\alpha }_n^{i}G\left(y_n{y}_n^i{V}_n^i\right)\mathrm{and}{\displaystyle \sum _{i=1}^N}{\alpha }_n^i:=1.$$

  

  This is graphically depicted in Figure 31 for y_n ∈

  

  . The interested reader is referred to [4] for more details.
• Substituting (17) into (16) yields the following p(y_n ∈ Y_n |x_n ) if t_n = t_ka : $$p\left(y_n\in Y_n|x_n\right)\approx {\displaystyle \sum _{i=1}^N}{\alpha }_n^i{\int }_{-\infty }^{\infty }G\left(y_n{\mathit{Cx}}_{n}V\right)G\left(y_n{y}_n^i{V}_n^i\right){\mathit{dy}}_n\mathrm{.}$$

  
• Proposition 1. [12,14] Let there exist two Gaussians of random vectors x ∈

  

  and m ∈

  

  , with Γ ∈

  

  ; G(m,Γx,M) and G(x,u,U). Then they satisfy: $${\int }_{-\infty }^{\infty }G\left(xuU\right)G\left(m\mathit{\Gamma x}M\right)\mathit{dx}=G\left(\mathit{\Gamma u},m,{\mathit{\Gamma U\Gamma }}^T+M\right),$$

  

  $$G\left(xuU\right)G\left(m\mathit{\Gamma x}M\right)=G\left(xdD\right)G\left(m,\mathit{\Gamma u},{\mathit{\Gamma U\Gamma }}^T+M\right),{}_{}\mathit{with\; D}:={\left(U^{-1}+{\mathit{\Gamma }}^T{M}^{-1}\mathrm{\Gamma }\right)}^{-1}{}_{}\mathit{and\; d}:={\mathit{DU}}^{-1}u+{\mathit{D\Gamma }}^T{M}^{-1}m\mathrm{.}$$

  
• Applying Proposition 1 ((19) to be precise) and G(x,y,Z) = G(y,x,Z) on (18) yields: $$p\left(y_n\in Y_n|x_n\right)\approx {\displaystyle \sum _{i=1}^N}{\alpha }_n^{i}G\left(y_n^i,{\mathit{Cx}}_n,V,+V_n^i\right),\mathrm{if}{t}_n=t_{k_a}\mathrm{.}$$

  
• In conclusion we can state that the unified expression of the PDF p(y_n ∈ Y_n |x_n ), at both t_n = t_ke and t_n = t_ka , for any event sampling method results in: $$p\left(y_n\in Y_n|x_n\right)\approx {\displaystyle \sum _{i=1}^N}{\alpha }_n^{i}G\left(y_n^i{\mathit{Cx}}_n{R}_n^i\right)\mathrm{with}{R}_n^i:=V+V_n^i\mathrm{.}$$

  
• If t_n = t_ke the variables of (22) are: N =1, ${\mathrm{<figure-callout\; id="1"\; label="\alpha \; n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">\alpha </figure-callout>}}_n^{}=1,{\mathrm{<figure-callout\; id="1"\; label="y\; n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">y</figure-callout>}}_n^{}=y_{k_e}$

  

  and ${\mathrm{<figure-callout\; id="1"\; label="V\; n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">V</figure-callout>}}_n^{}=\underline{0}\mathrm{.}$

  

  If t_n = t_ka the variables depend on ${\mathrm{\Lambda }}_{{\mathit{H}}_{{\mathit{k}}_{\mathit{e}}}}\left({\mathit{y}}_{\mathit{n}}\right)$

  

  and its approximation. As an example these variables are calculated for the method "Send-on-Delta" with y ∈

  

  .
• Example 1. In "Send-on-Delta", for certain N, the approximation of ${\mathrm{\Lambda }}_{{\mathit{H}}_{{\mathit{k}}_{\mathit{e}}}}\left({\mathit{y}}_{\mathit{n}}\right),$

  

  as presented in (17), is obtained with i ∈ {1, 2, ... ,N} and; $$\begin{array}{ll}{y}_n^i& =y_{k_e-1}-\left(\frac{N-2\left(i-1\right)-1}{2N}\right)2\mathrm{\Delta \; and}\\ {\alpha }_n^i& =1/N,V_n^i={\left(\frac{2\mathit{<figure-callout\; id="2"\; label="\Delta \; N"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">\Delta </figure-callout>}}{N}\right)}^2\left(0.25-0.05e^{-\frac{4\left(N-1\right)}{15}}-0.08e^{-\frac{4\left(N-1\right)}{150}}\right),\forall i\mathrm{.}\end{array}$$

  
• With the result of (22), p(x_n |y _0:n ∈ Y _0:n ) can also be expressed as a sum of N Gaussians.
5.2 Step 2: state estimation
• First the PDF p(x_n |y _0:n-1 ∈ Y _0:n-1) of (14b) is calculated. From the EBSE we have p(x _n-1|y_0:n-1 ∈ Y_0:n-1) := G(x _n-1,x _n-1|n-1, P _n-1, _n-1) and from (7a) with τ _n := t_n - t _n-1 we have $p\left(x_n|x_{n-1}\right):=G\left(x_n,A_{{\tau }_n}{x}_{n-1},B_{{\tau }_n}Q{B}_{{\tau }_n}^T\right)\mathrm{.}$

  

  Therefore using (19) in (14b) yields: $$p\left(x_n|y_{0:n-1}\in Y_{0:n-1}\right)=G\left(x_n{x}_{n|n-1}{P}_{n,n-1}\right)\mathrm{with}{x}_{n|n-1}:=A_{{\tau }_n}{x}_{n-1|n-1}\mathrm{and}{P}_{n|n-1}:=A_{{\tau }_n}{P}_{n-1\left|n-1\right|}{A}_{{\tau }_n}^T+B_{{\tau }_n}Q{B}_{{\tau }_n}^T\mathrm{.}$$

  
• Next p(x_n |y _0:n ∈ Y _0:n ), defined in (13), is calculated after multiplying (22) and (24): $$p\left(x_n|y_{n-1}\in Y_{0:n-1}\right)p\left(y_n\in Y_n\right)\approx {\displaystyle \sum _{i=1}^N}{\alpha }_n^{i}G\left(x_n{x}_{n|n-1}{P}_{n|n-1}\right)G\left(y_n^i{\mathit{Cx}}_n{R}_n^i\right)\mathrm{.}$$

  
• Equation (25) is explicitly solved by applying Proposition 1: $$p\left(x_n|y_{0:n-1}\in Y_{0:n-1}\right)p\left(y_n\in Y_n\right)\approx {\displaystyle \sum _{i=1}^N}{\alpha }_n^i{\beta }_n^{i}G\left(x_n{x}_n^i{P}_n^i\right)$$

  

  $$x_n^i:=P_n^i\left(P_{n|n-1}^{-1}{x}_{n|n-1}+C^T{\left(R_n^i\right)}^{-1}{y}_n^i\right),P_n^i:={\left(P_{n|n-1}^{-1}+C^T{\left(R_n^i\right)}^{-1}C\right)}^{-1}\mathrm{and}{\beta }_n^i:=G\left(y_n^i,{\mathit{Cx}}_{n|n-1},{\mathit{CP}}_{n|n-1}{C}^T+R_n^i\right)\mathrm{.}$$

  
• The expression of p(x_n |y _0:n ∈ Y _0:n ) as a sum of N Gaussians is the result of the following substitutions: (26) into (13), (26) into (14c) to obtain p(y_n ∈ Y_n | y _0:n-1 ∈ Y _0:n-1) and the latter into (13) again. This yields $$p\left(x_n|y_{0:n}\in Y_{0:n}\right)\approx {\displaystyle \sum _{i=1}^N}\frac{{\alpha }_n^i{\beta }_n^i}{\sum _{i=1}^N{\alpha }_n^i{\beta }_n^i}G\left(x_n{x}_n^i{P}_n^i\right)\mathrm{.}$$

  
• The third step is to approximate (27) as a single Gaussian to retrieve a computationally tractable algorithm. For if both p(x _n-1|y _0:n-1 ∈ Y _0:n-1) and p(y_n ∈ Y_n |x_n ) are approximated using N Gaussians, the estimate of x_n in (27) is described with M_n Gaussians. The value of M_n equals M _n-1 N, meaning that M_n increases after each sample instant and with it also the processing demand of the EBSE increases.
5.3 Step 3: state approximation
• p(x_n |y _0:n ∈ Y _0:n ) of (27) is approximated as a single Gaussian with an equal expectation and covariance matrix, i.e.: $$p\left(x_n|y_{0:n}\in Y_{0:n}\right)\approx G\left(x_n{x}_{n|n}{P}_{n|n}\right)$$

  

  $$x_{n|n}:={\displaystyle \sum _{i=1}^N}\frac{{\alpha }_n^i{\beta }_n^i{x}_n^i}{\sum _{i=1}^N{\alpha }_n^i{\beta }_n^i},P_{n|n}:={\displaystyle \sum _{i=1}^N}\frac{{\alpha }_n^i{\beta }_n^i}{\sum _{i=1}^N{\alpha }_n^i{\beta }_n^i}\left(P_n^i+\left(x_{n|n}-x_n^i\right){\left(x_{n|n}-x_n^i\right)}^T\right)\mathrm{.}$$

  
• The expectation an covariance of (27), equal to x _n|n and P _n|n of (28), can be derived from the corresponding definitions. Notice that because the designed EBSE is based on the equations of the Kalman filter, the condition of computational tractability is met.
5.4 Asymptotic analysis of the error-covariance matrix
• In this section we investigate the asymptotic analysis of the error-covariance matrix of the developed EBSE. By this we mean that we analyze lim_n→∞ P _n|n , which for convenience is denoted as P _∞. Note that for the classical Kalman filter (KF) [12] such an analysis is already available. However, for any other type of estimator asymptotic analysis remains a very challenging problem, which is why in most cases it is not even considered.
• Let us first recall the result on the asymptotic analysis of the Kalman filter. If x(t) of (7) is estimated, directly from y(t), with the KF at synchronous sampling times t_n := n · t_s , then P _n|n is updated as follows: $$P_{n|n}={\left({\left(A_{t_s}{P}_{n-1|n-1}{A}_{t_s}^T+B_{t_s}Q{B}_{t_s}^T\right)}^{-1}+C^T{V}^{-1}C\right)}^{-1}\mathrm{.}$$

  
• In [18,19] it is proven that if the eigenvalues of A_ts are within the unit circle and (A_ts ,C) is observable, then P _∞ = P_K . The matrix P_K equals the solution of: $$P_K={\left({\left(A_{t_s}{P}_K{A}_{t_s}^T+B_{t_s}Q{B}_{t_s}^T\right)}^{-1}+V^{-1}\right)}^{-1}\mathrm{.}$$

  

  For the EBSE however, we cannot prove that P _∞ equals a constant matrix. Instead we will prove that all the eigenvalues of P _∞ are bounded, i.e. that λ_max(P _∞) < ∞. As described in Section 4 this is a valid indication of an estimator's performance.
• The main result of this section is obtained under the standing assumption that ${\mathrm{\Lambda }}_{{\mathit{H}}_{{\mathit{k}}_{\mathit{e}}}}$

  

  is approximated using a single Gaussian. Note that the result then also applies to the estimator presented in [15], as a particular case. We assume that the eigenvalues of the A_τn-matrix are within the unit-circle and (A _τn ,C) is an observable pair. The following technical Lemmas will be of use.
• Lemma 1. Given the process model (7) and covariance matrices P > 0 and Q > 0, then for any 0 < τ₁ ≤ τ₂ we have that $A_{{\mathit{<figure-callout\; id="1"\; label="\tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">\tau </figure-callout>}}_1}P{A}_{{\mathit{\tau }}_1}^T\le A_{{\mathit{<figure-callout\; id="2"\; label="\tau "\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">\tau </figure-callout>}}_2}P{A}_{{\mathit{\tau }}_2}^T$
  
  and ${\mathrm{<figure-callout\; id="1"\; label="B\; \tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">B</figure-callout>}}_{{\mathit{\tau }}_{}}Q{B}_{{\mathit{\tau }}_1}^T\le {\mathrm{<figure-callout\; id="2"\; label="B\; \tau "\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">B</figure-callout>}}_{{\mathit{\tau }}_{}}Q{B}_{{\mathit{\tau }}_2}^T\mathrm{.}$
  
  See the Appendix for the proof.
• Lemma 2. Let any square matrices V₁ ≤ V₂ and W ₁ ≤ W₂ with V ₁ ≥ 0 and W ₁ ≥ 0 be given. Suppose that the matrices U ₁ and U ₂ are defined as $U_1:={\left(V_1^{-1}+C^T{W}_1^{-1}C\right)}^{-1}$
  
  and $U_2:={\left(V_2^{-1}+C^T{W}_2^{-1}C\right)}^{-1},$
  
  for any C of suitable size. Then it holds that U ₁ ≤ U ₂ .
• Proof. From [20] we have that $V_1^{-1}\ge V_2^{-1}$

  

  and $C^T{W}_1^{-1}C\ge C^T{W}_2^{-1}C\mathrm{.}$

  

  Hence, it follows that $V_1^{-1}+C^T{W}_1^{-1}C\ge V_2^{-1}+C^T{W}_2^{-1}C,$

  

  which yields $U_1^{-1}\ge U_2^{-1}\mathrm{.}$

  

  Thus, U ₁ ≤ U ₂, which concludes the proof.
• Next, recall that H_ke (y_n ) is assumed to be a bounded set. Therefore, it is reasonable to further assume that ${\mathrm{\Lambda }}_{{\mathit{H}}_{{\mathit{k}}_{\mathit{e}}}}$

  

  can be approximated using the formula (17), for N = 1, and that there exists a constant matrix V such that ${\mathrm{<figure-callout\; id="1"\; label="V\; n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">V</figure-callout>}}_n^{}\le \stackrel{}{V}$

  

  for all n.
• Theorem 1. Suppose that the EBSE, as presented in Section 5, approximates ${\mathrm{\Lambda }}_{{\mathit{H}}_{{\mathit{k}}_{\mathit{e}}}}$

  

  according to (17) with N = 1 and the above assumptions hold. Then λ _max (P _∞) < λ _max (P̃_K ), where P_K is equal to the solution of ${\stackrel{}{P}}_K={\left({\left(A_{t_s}{\stackrel{}{P}}_K{A}_{t_s}^T+B_{{\tau }_s}Q{B}_{{\tau }_s}^T\right)}^{-1}+{\left(V+\stackrel{}{V}\right)}^{-1}\right)}^{-1}\mathrm{.}$

  
• See the Appendix for the proof.
6 Illustrative example
• In this section we illustrate the effectiveness of the developed EBSE in terms of state-estimation error, sampling efficiency and computational tractability. The case study is a 1D object-tracking system. The states x(t) of the object are position and speed while the measurement vector y(t) is position. The process-noise w(t) represents the object's acceleration. Then given a maximum acceleration of 0.5[mls ²] its corresponding Q, according to [21], equals 0.02. Therefore the model as presented in (7) yields $A=\left(\begin{array}{cc}1& \mathrm{<figure-callout\; id="0"\; label="\tau "\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">\tau </figure-callout>}\\ 0& 1\end{array}\right),$

  

  $B=\left(\begin{array}{c}\frac{}{\mathrm{2}}\\ \mathrm{\tau }\end{array}\right),$

  

  C = (10) and D = 0, which is in fact a discrete-time double integrator. The acceleration in time is shown in Figure 32 together with the object's position and speed. The sampling time is t_s = 0.1 and the measurement-noise covariance is V = 0.1 · 10^-3.
• Three different estimators are tested. The first two estimators are the EBSE and the asynchronous Kalman filter (AKF) of [13]. For simplicity, in both estimators we used the "Send-on-Delta" method with Δ = 0.1 [nt]. For the EBSE we approximated ${\mathrm{\Lambda }}_{{\mathit{H}}_{{\mathit{k}}_{\mathit{e}}}}\left({\mathit{y}}_{\mathit{n}}\right)$

  

  using (23) with N = 5. The AKF estimates the states only at the event instants t_ke . The states at t_ka are calculated by applying the prediction-step of (14b). The third estimator is based on the quantized Kalman filter (QKF) introduced in [21] that uses synchronous time sampling af y_ka . The QKF can deal with quantized data, which also results in less data transfer, and therefore can be considered as an alternative to EBSE. In the QKF y _ka is the quantized version of y_ka with quantization level 0.1, which corresponds to the "Send-on-Delta" method. Hence, a comparison can be made.
• In Figure 33 and Figure 34 the state estimation-error of the three estimators is plotted. They show that the QKF estimates the position of the object with the least error. However, its error in speed is worse compared to the EBSE. Further, the plot of the AKF clearly shows that prediction of the state-estimates gives a significant growth in estimation-error when the time between the event sampling-instants increases (t > 4).
• Beside estimation error, sampling efficiency η is also important due to the increased interest in WSNs. For these systems communication is expensive and one aims to have the least data transfer. We define η ∈

  

  as $\eta :=\frac{{\left(x_i-x_{i|i}\right)}^T\left(x_i-x_{i|i}\right)}{{\left(x_i-x_{i|i-1}\right)}^T\left(x_i-x_{i|i-1}\right)},$

  

  which is a measure of the change in the estimation-error after the measurement update with either z_ke or y _ka was done. Notice that if η < 1 the estimation error decreased after an update, if η > 1 the error increased and if η = 1 the error remained the same. For the EBSE i = k_e with i - 1 equal to k_e - 1 or k_a - 1. For the AKF i = k_e with i - 1 = k_e - 1. For the QKF i = k_a and i - 1 = k_a - 1. Figure 35 shows that for the EBSE η < 1 at all instants n. The AKF has one instant, t = 3.4, at which η > 1. In case of the QKF the error sometimes decreases but it can also increase considerably after an update. Also notice that η of the QKF converges to 1. Meaning that for t > 5.5 the estimation error does not change after an update and new samples are mostly used to bound λ _i (P _{ka |ka} ). The EBSE has the same property, although for this method the last sample was received at t = 4.9.
• The last aspect on which the three estimators are compared is the total amount of processing time which was needed to estimate all state-vectors. For the EB SE, both x_ke and x_ka were estimated and it took 0.094 seconds. The AKF estimated x_ke and predicted x_ka in a total time of 0.016 seconds and the QKF estimated x_ka and its total processing time equaled 0.022 seconds. This means that although the EBSE results in the most processing time, it is computationally comparable to the AKF and QKF, while it provides an estimation-error similar to the QKF, but with significantly less data transmission. As such, it is most suited for usage in networks in general and WSNs in particular.
7 Conclusions
• In this paper a general event-based state-estimator was presented. The distinguishing feature of the proposed EBSE is that estimation of the states is performed at two different type of time instants, i.e. at event instants t_ke , when measurement data is used for update, and at synchronous time sampling t_ka , when no measurement is received, but an update is performed based on the knowledge that the monitored variable lies within a set used to define the event. As a result, it could be proven that, under certain assumptions, for the error-covariance matrix of the EBSE it holds that λ _max (P _∞) < ∞, even in the situation when no new observation z_ke is received anymore. Its effectiveness for usage in WSNs has been demonstrated on an application example.
References
• 1. L F. Akyildiz, W. Su, Y Sankarasubramaniam, and E. Cayirci, "Wireless Sensor Networks: a survey," Elsevier, Computer Networks, vol. 38, pp. 393-422, 2002.
• 2. E. Johannesson, T. Henningsson, and A. Cervin, "Sporadic control of first-order linear stochastic systems," in Hybrid Systems: Computation and Control, ser. Lecture Notes in Computer Science, vol. 4416. Pisa, Italy: Springer Verlag, 2007, pp. 301-314.
• 3. W. P. M. H. Heemels, J. H. Sandee, and P. P. J. van den Bosch, "Analysis of event-driven controllers for linear systems," International Journal of Control, vol. 81, no. 4, 2008.
• 4. H. W. Sorenson and D. L. Alspach, "Recursive Bayesian estimation using Gaussian sums," Automatica, vol. 7, pp. 465-479, 1971.
• 5. J. H. Kotecha and P. M. Djurié, "Gaussian sum particle filtering," IEEE Transaction Signal Processing, vol. 51, no. 10, pp. 2602-2612,2003.
• 6. N. L. Johnson, S. Kotz, and A. W. Kemp, Univariate discrete distributions, John Wiley and Sons, 1992.
• 7. L. Aggoun and R. Elliot, Measure Theory and Filtering. Cambridge University Press, 2004.
• 8. H. L. Lebesque, "Integrale, longueur, aire," Ph.D. dissertation, University of Nancy, 1902.
• 9. K. J. Åström and B. M. Bernhardsson, "Comparison of Riemann and Lebesque sampling for first order stochastic systems," in 41st IEEE Conf. on Dec. and Contr., Las Vegas, USA, 2002.
• 10. M. Miskowicz, "Send-on-delta concept: an event-based data-reporting strategy," Sensors, vol. 6, pp. 49-63,2006.
• 11. ----, "Asymptotic Effectiveness of the Event-Based Sampling according to the Integral Criterion," Sensors, vol. 7, pp. 16-37, 2007.
• 12. R. E. Kalman, "A new approach to linear filtering and prediction problems," Transaction of the ASME Journal of Basic Engineering, vol. 82, no. D, pp. 35-42, 1960.
• 13. M. Mallick, S. Coraluppi, and C. Carthel, "Advances in Asynchronous and Decentralized Estimation," in In Proceeding of the 2001 Aerospace Conference, Big Sky, MT, USA, 2001.
• 14. B. Ristic, S. Arulampalam, and N. Gordon, Beyond the Kalman filter: Particle filter for tracking applications. Boston, Massachusetts: Artech House, 2004.
• 15. V. H. Nguyen and Y S. Suh, "Improving estimation performance in Networked Control Systems applying the Send-on-delta transmission method," Sensors, vol. 7, pp. 2128-2138, 2007.
• 16. K. V. Mardia, J. T. Kent, and J. M. Bibby, Mutlivariate analysis. Academic press, London, 1979.
• 17. D. C. Montgomery and G. C. Runger, Applied Statistics and Probability for Engineers. John Wiley and Sons, 2007.
• 18. M. L.J. Hautus, "Controllability and observability conditions of linear autonomous systems," in Indagationes Mathemathicae, vol. 32, 1972, pp. 448-455.
• 19. A. V Balakrishnan, Kalman Filtering Theory. Optimization Software, New York, 1987.
• 20. D. S. Bernstein, Matrix Mathematics. Princeton University Press, 2005.
• 21. R. E. Curry, Estimation and control with quantized measurements. Boston, Massachusetts: MIT Press, 1970.
• 22. G. F. Franklin, J. D. Powel, and A. Emami-Naeini, Feedback control of dynamic systems. Addison-Wesley 1995.
A Proof of Lemma 1
• Suppose that A ∈

  

  and B ∈

  

  are defined as the state-space matrices for the time-continuous counterpart of (7). Then it is known [22] that for any sampling period τ > 0, A _τ and B _τ of (7a) are obtained from their corresponding continuous-time matrices A and B as follows: $$A_{\tau }:=e^{\mathit{A\tau }}:={\displaystyle \sum _{i=0}^{\infty }}{A}^i\frac{{\tau }^i}{i!}\mathrm{and}{B}_{\tau }:={\int }_0^{\tau }{e}^{\mathit{A\eta }}\mathit{d\eta B}:={\displaystyle \sum _{i=0}^{\infty }}{A}^{i}B\frac{{\tau }^{i+1}}{\left(i+1\right)!}\mathrm{.}$$

  
• Using (31) one obtains: $$A_{{\mathrm{<figure-callout\; id="2"\; label="\tau "\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">\tau </figure-callout>}}_2}{\mathit{PA}}_{{\tau }_2}^T-{\mathrm{<figure-callout\; id="1"\; label="A\; \tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">A</figure-callout>}}_{{\tau }_{}}{\mathit{PA}}_{{\tau }_1}^T=\left({\displaystyle \sum _{i=0}^{\infty }}{A}^i\frac{{\tau }_2^i}{i!}\right)P\left({\displaystyle \sum _{j=0}^{\infty }}{\left(A^T\right)}^j\frac{{\tau }_2^i}{j!}\right)-\left({\displaystyle \sum _{i=0}^{\infty }}{A}^i\frac{{\tau }_1^i}{i!}\right)P\left({\displaystyle \sum _{j=0}^{\infty }}{\left(A^T\right)}^j\frac{{\tau }_1^j}{j!}\right)={\displaystyle \sum _{i=0}^{\infty }}{\displaystyle \sum _{j=0}^{\infty }}{A}^{i}P{\left(A^T\right)}^j\left(\frac{{\tau }_2^i}{i!}\frac{{\tau }_2^j}{j!}\right)-{\displaystyle \sum _{i=0}^{\infty }}{\displaystyle \sum _{j=0}^{\infty }}{A}^{i}P{\left(A^T\right)}^j\left(\frac{{\tau }_1^i}{i!}\frac{{\tau }_1^j}{j!}\right)={\displaystyle \sum _{i=0}^{\infty }}{\displaystyle \sum _{j=0}^{\infty }}{A}^{i}P{\left(A^T\right)}^j\left(\frac{{\tau }_2^{i+j}-{\tau }_1^{i+j}}{i!j!}\right)$$

  
• As for any τ > 0 the series e ^Aτ converges [22], then $A_{{\mathit{<figure-callout\; id="2"\; label="\tau "\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">\tau </figure-callout>}}_2}P{A}_{{\mathit{\tau }}_2}^T-{\mathrm{<figure-callout\; id="1"\; label="A\; \tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">A</figure-callout>}}_{{\mathit{\tau }}_{}}P{A}_{{\mathit{\tau }}_1}^T$

  

  also converges. Then, since 0 < τ₁ ≤ τ₂ and P > 0, for any fixed i, j, we have $A^{i}P{\left(A^T\right)}^j\left(\frac{{\mathrm{\tau }}_2^{i+j}-{\mathrm{\tau }}_1^{i+j}}{i\left|j\right|}\right)\ge 0$

  

  for any matrix A and thus, it follows that $A_{{\mathit{<figure-callout\; id="2"\; label="\tau "\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">\tau </figure-callout>}}_2}P{A}_{{\mathit{\tau }}_2}^T\ge {\mathrm{<figure-callout\; id="1"\; label="A\; \tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">A</figure-callout>}}_{{\mathit{\tau }}_{}}P{A}_{{\mathit{\tau }}_1}^T\mathrm{.}$

  

  The same reasoning can be used to prove that ${\mathrm{<figure-callout\; id="1"\; label="B\; \tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">B</figure-callout>}}_{{\mathit{\tau }}_{}}Q{B}_{{\mathit{\tau }}_1}^T\le {\mathrm{<figure-callout\; id="2"\; label="B\; \tau "\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">B</figure-callout>}}_{{\mathit{\tau }}_{}}Q{B}_{{\mathit{\tau }}_2}^T,$

  
B Proof of Theorem 1
• Under the hypothesis, for the proposed EBSE, P _n|n of (28), with τ _n := t_n - t _n-1 and $R_n:=V+{\mathrm{<figure-callout\; id="1"\; label="V\; n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">V</figure-callout>}}_n^{},$

  

  becomes: $$P_{n|n}={\left({\left(A_{{\tau }_n}{P}_{n-1|n-1}{A}_{{\tau }_n}^T+B_{{\tau }_n}Q{B}_{{\tau }_n}^T\right)}^{-1}+C^T{R}_n^{-1}C\right)}^{-1}\mathrm{.}$$

  
• The upper bound on λ _max (P _∞) is proven by induction, considering the asymptotic behavior of a KF that runs in parallel with the EBSE, as follows. The EBSE calculates $P_{n|\mathrm{<figure-callout\; id="1"\; label="n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(1\right)}$

  

  as (32) and the KF calculates $P_{n|\mathrm{<figure-callout\; id="2"\; label="n"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(2\right)}$

  

  as (29) in which V is replaced with R := V + V . Notice that for these estimators we have that τ _n ≤ t_s and R_n ≤ R, for all n. Let the EBSE and the KF start with the same initial covariance matrix P ₀.
• The first step of induction is to prove that ${\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(1\right)}\le {\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(2\right)}\mathrm{.}$

  

  From the definition of ${\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(1\right)}$

  

  in (32) and ${\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(2\right)}$

  

  in (29) we have that ${\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(1\right)}={\left({\left(A_{{\mathrm{<figure-callout\; id="1"\; label="\tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">\tau </figure-callout>}}_1}{\mathrm{<figure-callout\; id="0"\; label="P"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">P</figure-callout>}}_0{A}_{{\tau }_1}^T+{\mathrm{<figure-callout\; id="1"\; label="B\; \tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">B</figure-callout>}}_{{\tau }_{}}Q{B}_{{\tau }_1}^T\right)}^{-1}+C^T{R}_1^{-1}C\right)}^{-1}$

  

  and ${\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(2\right)}:={\left({\left(A_{t_s}{\mathrm{<figure-callout\; id="0"\; label="P"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">P</figure-callout>}}_0{A}_{t_s}^T+B_{t_s}Q{B}_{t_s}^T\right)}^{-1}+C^T{R}^{-1}C\right)}^{-1}\mathrm{.}$

  
• Suppose we define $V_1:=A_{{\mathrm{<figure-callout\; id="1"\; label="\tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">\tau </figure-callout>}}_1}{\mathrm{<figure-callout\; id="0"\; label="P"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">P</figure-callout>}}_0{A}_{{\tau }_1}^T+{\mathrm{<figure-callout\; id="1"\; label="B\; \tau "\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">B</figure-callout>}}_{{\tau }_{}}Q{B}_{{\tau }_1}^T,$

  

  $V_2:=A_{t_s}{\mathrm{<figure-callout\; id="0"\; label="P"\; filenames="imgf0017.png,imgf0018.png"\; state="\{\{state\}\}">P</figure-callout>}}_0{A}_{t_s}^T+B_{t_s}Q{B}_{t_s}^T,$

  

  W ₁ := R ₁ and W ₂ := R. Then W ₁ ≤ W ₂ and from Lemma 1 it follows that V ₁ ≤ V ₂. Therefore applying Lemma 2, with $U_1:={\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(1\right)}$

  

  and $U_2:={\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(2\right)},$

  

  yields ${\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(1\right)}\le {\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">P</figure-callout>}}_{1|1}^{\left(2\right)}\mathrm{.}$

  
• The second and last step of induction is to show that if $P_{n-1|n-1}^{\left(1\right)}\le P_{n-1|n-1}^{\left(2\right)},$

  

  then $P_{n|\mathrm{<figure-callout\; id="1"\; label="n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(1\right)}\le P_{n|\mathrm{<figure-callout\; id="2"\; label="n"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(2\right)}\mathrm{.}$

  

  Let $V_1:=A_{{\tau }_n}{P}_{n-1|n-1}^{\left(1\right)}{A}_{{\tau }_n}^T+B_{{\tau }_n}Q{B}_{{\tau }_n}^T,$

  

  $V_2:=A_{t_s}{P}_{n-1|n-1}^{\left(2\right)}{A}_{t_s}^T+B_{t_s}Q{B}_{t_s}^T,$

  

  and let W ₁ := R_n , W ₂ := R. Notice that this yields W ₁ ≤ W ₂ . The second condition of Lemma 2, i.e. V ₁ ≤ V ₂ also holds by applying Lemma 1, i.e. $$A_{{\tau }_n}{P}_{n-1|n-1}^{\left(1\right)}{A}_{{\tau }_n}^T+B_{{\tau }_n}Q{B}_{{\tau }_n}^T\le A_{t_s}{P}_{n-1|n-1}^{\left(1\right)}{A}_{t_s}^T+B_{t_s}Q{B}_{t_s}^T\le A_{t_s}{P}_{n-1|n-1}^{\left(2\right)}{A}_{t_s}^T+B_{t_s}Q{B}_{t_s}^T\mathrm{.}$$

  

  Hence, applying Lemma 2, with $U_1:=P_{n|\mathrm{<figure-callout\; id="1"\; label="n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(1\right)}$

  

  and $U_2:=P_{n|\mathrm{<figure-callout\; id="2"\; label="n"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(2\right)}$

  

  yields $P_{n|\mathrm{<figure-callout\; id="1"\; label="n"\; filenames="imgf0008.png,imgf0011.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(1\right)}\le P_{n|\mathrm{<figure-callout\; id="2"\; label="n"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(2\right)}\mathrm{.}$

  
• This proves that $P_{\infty }^{\left(1\right)}\le P_{\infty }^{\left(2\right)},$

  

  which yields (see e.g., [20]) ${\lambda }_{\mathit{max}}\left(P_{\infty }^{\left(1\right)}\right)\le {\lambda }_{\mathit{max}}\left(P_{\infty }^{\left(2\right)}\right)\mathrm{.}$

  

  As $P_{n|\mathrm{<figure-callout\; id="2"\; label="n"\; filenames="imgf0019.png,imgf0021.png"\; state="\{\{state\}\}">n</figure-callout>}}^{\left(2\right)}$

  

  was calculated with the KF it follows from (30) that $P_{\infty }^{\left(2\right)}={\stackrel{}{P}}_K,$

  

  with P̃_K as the solution of ${\stackrel{}{P}}_K={\left({\left(A_{t_s}{\stackrel{}{P}}_K{A}_{t_s}^T+B_{t_s}Q{B}_{t_s}^T\right)}^{-1}+R^{-1}\right)}^{-1},$

  

  which completes the proof.

Claims (12)

1. Traffic information unit (MD1, MD2, MD3,...) associated with a traffic infrastructure comprising

   - a facility (MI) for tracking vehicle state information of individual vehicles present at the traffic infrastructure,

   - a facility (T) for transmitting said vehicle state information to a vehicle (70B, 70E) at the traffic infrastructure.
2. Traffic information unit (MD1, MD2, MD3,...) further comprising

   - a sensor system comprising a plurality of sensor nodes (10) for sensing vehicles (70A, ...,70E) arranged in the vicinity of a traffic infrastructure (80) for carrying vehicles,

   - communication means (16) coupled to the sensors, wherein the facility (MI) is a message interpreter that uses information (D) communicated by the sensor nodes.
3. Traffic information unit (MD1, MD2, MD3) according to claim 2, wherein the sensor nodes (10) provide a message (D) indicative for an occupancy status of a detection area of a traffic infrastructure monitored by the sensor nodes, the message interpreter (MI) further comprising:

   - a vehicle database facility (32, 34) comprising state information of vehicles present at the traffic infrastructure, the vehicle state information including at least one of a vehicle position, a vehicle speed, a vehicle orientation,

   - an association facility (40) for associating the messages (D) provided by the sensor elements (10) with the vehicle state information present in the vehicle data base facility,

   - a state updating facility (50) for updating the vehicle state information on the basis of the messages associated therewith.
4. Traffic information unit according to claim 2, wherein the sensor nodes (10) provide spatial occupancy information with a density higher than 1 m^-2.
5. Traffic information unit according to claim 2, wherein the sensor nodes (10) transmit data upon detection of an event.
6. Traffic information unit according to claim 2, wherein the sensor nodes (10) are embedded in the traffic infrastructure.
7. Traffic information system comprising at least a first and a second traffic information unit (MD1, MD2, MD3) according to one of the previous claims, the first and the second traffic information unit being associated with mutually neighboring sections of the traffic infrastructure and being arranged to mutually exchange vehicle state information.
8. A traffic information system according to claim 7, further comprising at least one client information module (CIM) for providing status information related to the infrastructure (80), the status comprising at least one of an occupation density and an average speed as a function of a position at the traffic infrastructure
9. A vehicle management system (C) for a target vehicle (70B, 70E) comprising a communication system (R) arranged for receiving vehicle state information relating to surrounding vehicles from a traffic information unit (MD1, MD2, MD3) according to one of the claims 1-6 or from a traffic information system according to claim 7 or 8, inputs (C1) for receiving vehicle state information from the target vehicle (70B, 70E) and a control system (C2) with outputs (C3) for providing control signals for controlling a state of the vehicle using the vehicle state information retrieved from the traffic information system.
10. A vehicle management system (C) according to claim 9, further comprising communication means (R1) for exchanging vehicle state information with surrounding vehicles and a selection facility (SL) for selecting one or more of vehicle state information obtained from the surrounding vehicles (VS2) and vehicle state information received (VS1) from the traffic information system as the vehicle state information (VS) to be used by the control system (C2).
11. A vehicle (70B, 70E) comprising a vehicle management system (C) according to claim 9 or 10.
12. Method of controlling a vehicle, comprising the steps of

    - observing vehicles from a fixed position,

    - communicating the observations,

    - tracking motion states of individual vehicles using the communicated observations

    - transmitting said information about said tracked states to a vehicle instrumented with a vehicle management system according to claim 9.

Images