CustusX  16.5
An IGT application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
cxTrackingPositionFilter.cpp
Go to the documentation of this file.
1 /*=========================================================================
2 This file is part of CustusX, an Image Guided Therapy Application.
3 
4 Copyright (c) 2008-2014, SINTEF Department of Medical Technology
5 All rights reserved.
6 
7 Redistribution and use in source and binary forms, with or without
8 modification, are permitted provided that the following conditions are met:
9 
10 1. Redistributions of source code must retain the above copyright notice,
11  this list of conditions and the following disclaimer.
12 
13 2. Redistributions in binary form must reproduce the above copyright notice,
14  this list of conditions and the following disclaimer in the documentation
15  and/or other materials provided with the distribution.
16 
17 3. Neither the name of the copyright holder nor the names of its contributors
18  may be used to endorse or promote products derived from this software
19  without specific prior written permission.
20 
21 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
22 AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
23 IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
24 DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
25 FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
26 DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
27 SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28 CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
29 OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
30 OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31 =========================================================================*/
33 #include "iir/Butterworth.h"
34 
35 namespace cx
36 {
37 
39 {
40  mCutOffFrequency = 3;
41  mResampleFrequency = 100;
42 
43  fx.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency); // Lag perker isteden
44  fx.reset ();
45  fy.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency);
46  fy.reset ();
47  fz.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency);
48  fz.reset ();
49 }
50 
52 {
53 
54  this->clearIfTimestampIsOlderThanHead(pos, timestamp);
55  this->clearIfJumpInTimestamps(pos, timestamp);
56 
57  if (mResampled.empty()){
58  mResampled[timestamp] = pos;
59  mHistory[timestamp] = pos;
60  return;
61  }
62 
63  this->interpolateAndFilterPositions(pos, timestamp);
64  mHistory[timestamp] = pos;
65 }
66 
68 {
69  if (mFiltered.size() > mResampleFrequency) //check if mFiltered contains enough positions for the filter to be stable
70  return mFiltered.rbegin()->second;
71  else if (!mHistory.empty())
72  return mHistory.rbegin()->second;
73  else
74  return Transform3D::Identity();
75 }
76 
77 void TrackingPositionFilter::clearIfTimestampIsOlderThanHead(Transform3D pos, double timestamp)
78 {
79  if (mResampled.empty())
80  return;
81 
82  if (timestamp < mResampled.rbegin()->first){ // clear history if old timestamps appear
83  mHistory.clear();
84  mResampled.clear();
85  mFiltered.clear();
86 
87  fx.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency); // Lag perker isteden
88  fx.reset ();
89  fy.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency);
90  fy.reset ();
91  fz.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency);
92  fz.reset ();
93  }
94 
95 }
96 
97 void TrackingPositionFilter::clearIfJumpInTimestamps(Transform3D pos, double timestamp)
98 {
99  if (mResampled.empty())
100  return;
101 
102  double timeStep = timestamp - mResampled.rbegin()->first;
103  if ( timeStep > 1000){ // clear history of resampled and filtered data if jump in timestamps of more than 1 second
104  mHistory.clear();
105  mResampled.clear();
106  mFiltered.clear();
107 
108  fx.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency); // Lag perker isteden
109  fx.reset ();
110  fy.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency);
111  fy.reset ();
112  fz.setup (mFilterOrder, mResampleFrequency, mCutOffFrequency);
113  fz.reset ();
114  }
115 }
116 
117 void TrackingPositionFilter::interpolateAndFilterPositions(Transform3D pos, double timestamp)
118 {
119  Transform3D previousPositionMatrix = mHistory.rbegin()->second;
120  double deltaT = timestamp - mHistory.rbegin()->first; //time from previous measured position to this position
121  int numberOfInterpolationPoints = floor( (timestamp - mResampled.rbegin()->first)/1000 * mResampleFrequency ); // interpolate from last resampled position to current measured position
122  Transform3D interpolatedPosition;
123  Transform3D filteredPosition;
124  for (int i=0; i < numberOfInterpolationPoints; i++)
125  {
126  double resampledTimestamp = mResampled.rbegin()->first + 1000/mResampleFrequency;
127  double deltaTpast = resampledTimestamp - mHistory.rbegin()->first;
128  double deltaTfuture = timestamp - resampledTimestamp;
129  interpolatedPosition = pos.matrix() * deltaTpast/deltaT + previousPositionMatrix.matrix() * deltaTfuture/deltaT; // linear interpolation between previous and current measured position
130  mResampled[resampledTimestamp] = interpolatedPosition;
131 
132  filteredPosition = interpolatedPosition;
133  filteredPosition(0,3) = fx.filter(interpolatedPosition(0,3));
134  filteredPosition(1,3) = fy.filter(interpolatedPosition(1,3));
135  filteredPosition(2,3) = fz.filter(interpolatedPosition(2,3));
136  mFiltered[resampledTimestamp] = filteredPosition;
137  }
138 
139 }
140 
141 } // namespace cx
142 
143 
Transform3D Transform3D
Transform3D is a representation of an affine 3D transform.
void addPosition(Transform3D pos, double timestamp)