| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  | """
 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  | GTSAM Copyright 2010-2020, Georgia Tech Research Corporation, | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  | Atlanta, Georgia 30332-0415 | 
					
						
							|  |  |  | All Rights Reserved | 
					
						
							|  |  |  | Authors: Frank Dellaert, et al. (see THANKS for the full author list) | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | See LICENSE for the license information | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | Track a moving object "Time of Arrival" measurements at 4 microphones. | 
					
						
							|  |  |  | Author: Frank Dellaert | 
					
						
							|  |  |  | """
 | 
					
						
							|  |  |  | # pylint: disable=invalid-name, no-name-in-module | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | from gtsam import (LevenbergMarquardtOptimizer, LevenbergMarquardtParams, | 
					
						
							| 
									
										
										
										
											2020-08-18 23:02:35 +08:00
										 |  |  |                    NonlinearFactorGraph, Point3, Values, noiseModel) | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  | from gtsam_unstable import Event, TimeOfArrival, TOAFactor | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | # units | 
					
						
							|  |  |  | MS = 1e-3 | 
					
						
							|  |  |  | CM = 1e-2 | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | # Instantiate functor with speed of sound value | 
					
						
							|  |  |  | TIME_OF_ARRIVAL = TimeOfArrival(330) | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | 
 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  | def define_microphones(): | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     """Create microphones.""" | 
					
						
							|  |  |  |     height = 0.5 | 
					
						
							|  |  |  |     microphones = [] | 
					
						
							|  |  |  |     microphones.append(Point3(0, 0, height)) | 
					
						
							|  |  |  |     microphones.append(Point3(403 * CM, 0, height)) | 
					
						
							|  |  |  |     microphones.append(Point3(403 * CM, 403 * CM, height)) | 
					
						
							|  |  |  |     microphones.append(Point3(0, 403 * CM, 2 * height)) | 
					
						
							|  |  |  |     return microphones | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | 
 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  | def create_trajectory(n): | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     """Create ground truth trajectory.""" | 
					
						
							|  |  |  |     trajectory = [] | 
					
						
							|  |  |  |     timeOfEvent = 10 | 
					
						
							|  |  |  |     # simulate emitting a sound every second while moving on straight line | 
					
						
							|  |  |  |     for key in range(n): | 
					
						
							|  |  |  |         trajectory.append( | 
					
						
							|  |  |  |             Event(timeOfEvent, 245 * CM + key * 1.0, 201.5 * CM, (212 - 45) * CM)) | 
					
						
							|  |  |  |         timeOfEvent += 1 | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  |     return trajectory | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | 
 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  | def simulate_one_toa(microphones, event): | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     """Simulate time-of-arrival measurements for a single event.""" | 
					
						
							|  |  |  |     return [TIME_OF_ARRIVAL.measure(event, microphones[i]) | 
					
						
							|  |  |  |             for i in range(len(microphones))] | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | 
 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  | def simulate_toa(microphones, trajectory): | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     """Simulate time-of-arrival measurements for an entire trajectory.""" | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  |     return [simulate_one_toa(microphones, event) | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |             for event in trajectory] | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | 
 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  | def create_graph(microphones, simulatedTOA): | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     """Create factor graph.""" | 
					
						
							|  |  |  |     graph = NonlinearFactorGraph() | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  |     # Create a noise model for the TOA error | 
					
						
							| 
									
										
										
										
											2020-08-18 23:02:35 +08:00
										 |  |  |     model = noiseModel.Isotropic.Sigma(1, 0.5 * MS) | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  | 
 | 
					
						
							|  |  |  |     K = len(microphones) | 
					
						
							|  |  |  |     key = 0 | 
					
						
							|  |  |  |     for toa in simulatedTOA: | 
					
						
							|  |  |  |         for i in range(K): | 
					
						
							|  |  |  |             factor = TOAFactor(key, microphones[i], toa[i], model) | 
					
						
							|  |  |  |             graph.push_back(factor) | 
					
						
							|  |  |  |         key += 1 | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  |     return graph | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | 
 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  | def create_initial_estimate(n): | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     """Create initial estimate for n events.""" | 
					
						
							|  |  |  |     initial = Values() | 
					
						
							|  |  |  |     zero = Event() | 
					
						
							|  |  |  |     for key in range(n): | 
					
						
							|  |  |  |         TOAFactor.InsertEvent(key, zero, initial) | 
					
						
							|  |  |  |     return initial | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | 
 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  | def toa_example(): | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     """Run example with 4 microphones and 5 events in a straight line.""" | 
					
						
							|  |  |  |     # Create microphones | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  |     microphones = define_microphones() | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     K = len(microphones) | 
					
						
							|  |  |  |     for i in range(K): | 
					
						
							|  |  |  |         print("mic {} = {}".format(i, microphones[i])) | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  |     # Create a ground truth trajectory | 
					
						
							|  |  |  |     n = 5 | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  |     groundTruth = create_trajectory(n) | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     for event in groundTruth: | 
					
						
							|  |  |  |         print(event) | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  |     # Simulate time-of-arrival measurements | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  |     simulatedTOA = simulate_toa(microphones, groundTruth) | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     for key in range(n): | 
					
						
							|  |  |  |         for i in range(K): | 
					
						
							|  |  |  |             print("z_{}{} = {} ms".format(key, i, simulatedTOA[key][i] / MS)) | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  |     # create factor graph | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  |     graph = create_graph(microphones, simulatedTOA) | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     print(graph.at(0)) | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  |     # Create initial estimate | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  |     initial_estimate = create_initial_estimate(n) | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     print(initial_estimate) | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  |     # Optimize using Levenberg-Marquardt optimization. | 
					
						
							|  |  |  |     params = LevenbergMarquardtParams() | 
					
						
							|  |  |  |     params.setAbsoluteErrorTol(1e-10) | 
					
						
							|  |  |  |     params.setVerbosityLM("SUMMARY") | 
					
						
							|  |  |  |     optimizer = LevenbergMarquardtOptimizer(graph, initial_estimate, params) | 
					
						
							|  |  |  |     result = optimizer.optimize() | 
					
						
							|  |  |  |     print("Final Result:\n", result) | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | 
 | 
					
						
							|  |  |  | if __name__ == '__main__': | 
					
						
							| 
									
										
										
										
											2020-03-25 20:59:58 +08:00
										 |  |  |     toa_example() | 
					
						
							| 
									
										
										
										
											2020-03-19 05:28:40 +08:00
										 |  |  |     print("Example complete") |