GestureRecognitionToolkit  Version: 1.0 Revision: 04-03-15
The Gesture Recognition Toolkit (GRT) is a cross-platform, open-source, c++ machine learning library for real-time gesture recognition.
DTW.cpp
1 /*
2  GRT MIT License
3  Copyright (c) <2012> <Nicholas Gillian, Media Lab, MIT>
4 
5  Permission is hereby granted, free of charge, to any person obtaining a copy of this software
6  and associated documentation files (the "Software"), to deal in the Software without restriction,
7  including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense,
8  and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so,
9  subject to the following conditions:
10 
11  The above copyright notice and this permission notice shall be included in all copies or substantial
12  portions of the Software.
13 
14  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT
15  LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
16  IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
17  WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
18  SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
19  */
20 
21 #include "DTW.h"
22 
23 using namespace std;
24 
25 namespace GRT{
26 
27 //Register the DTW module with the Classifier base class
28 RegisterClassifierModule< DTW > DTW::registerModule("DTW");
29 
30 DTW::DTW(bool useScaling,bool useNullRejection,double nullRejectionCoeff,UINT rejectionMode,bool constrainWarpingPath,double radius,bool offsetUsingFirstSample,bool useSmoothing,UINT smoothingFactor)
31 {
32  this->useScaling=useScaling;
33  this->useNullRejection = useNullRejection;
34  this->nullRejectionCoeff = nullRejectionCoeff;
35  this->rejectionMode = rejectionMode;
36  this->constrainWarpingPath = constrainWarpingPath;
37  this->radius = radius;
38  this->offsetUsingFirstSample = offsetUsingFirstSample;
39  this->useSmoothing = useSmoothing;
40  this->smoothingFactor = smoothingFactor;
41 
42  supportsNullRejection = true;
43  trained=false;
44  useZNormalisation=false;
45  constrainZNorm=false;
46  trimTrainingData = false;
47 
48  zNormConstrainThreshold=0.2;
49  trimThreshold = 0.1;
50  maximumTrimPercentage = 90;
51 
52  numTemplates=0;
53  distanceMethod=EUCLIDEAN_DIST;
54 
55  averageTemplateLength =0;
56 
57  classType = "DTW";
58  classifierType = classType;
59  classifierMode = TIMESERIES_CLASSIFIER_MODE;
60  debugLog.setProceedingText("[DEBUG NDDTW]");
61  errorLog.setProceedingText("[ERROR NDDTW]");
62  trainingLog.setProceedingText("[TRAINING NDDTW]");
63  warningLog.setProceedingText("[WARNING NDDTW]");
64 }
65 
66 DTW::DTW(const DTW &rhs){
67  *this = rhs;
68 }
69 
70 DTW::~DTW(void)
71 {
72 }
73 
74 DTW& DTW::operator=(const DTW &rhs){
75 
76  if( this != &rhs ){
77 
78  this->templatesBuffer = rhs.templatesBuffer;
79  this->distanceMatrices = rhs.distanceMatrices;
80  this->warpPaths = rhs.warpPaths;
81  this->continuousInputDataBuffer = rhs.continuousInputDataBuffer;
82  this->numTemplates = rhs.numTemplates;
83  this->rejectionMode = rhs.rejectionMode;
84  this->useSmoothing = rhs.useSmoothing;
85  this->useZNormalisation = rhs.useZNormalisation;
86  this->constrainZNorm = rhs.constrainZNorm;
87  this->constrainWarpingPath = rhs.constrainWarpingPath;
88  this->trimTrainingData = rhs.trimTrainingData;
89  this->zNormConstrainThreshold = rhs.zNormConstrainThreshold;
90  this->radius = rhs.radius;
91  this->offsetUsingFirstSample = rhs.offsetUsingFirstSample;
92  this->trimThreshold = rhs.trimThreshold;
93  this->maximumTrimPercentage = rhs.maximumTrimPercentage;
94  this->smoothingFactor = rhs.smoothingFactor;
95  this->distanceMethod = rhs.distanceMethod;
96  this->rejectionMode = rhs.rejectionMode;
97  this->averageTemplateLength = rhs.averageTemplateLength;
98 
99  //Copy the classifier variables
100  copyBaseVariables( (Classifier*)&rhs );
101  }
102  return *this;
103 }
104 
105 bool DTW::deepCopyFrom(const Classifier *classifier){
106 
107  if( classifier == NULL ) return false;
108 
109  if( this->getClassifierType() == classifier->getClassifierType() ){
110 
111  DTW *ptr = (DTW*)classifier;
112  this->templatesBuffer = ptr->templatesBuffer;
113  this->distanceMatrices = ptr->distanceMatrices;
114  this->warpPaths = ptr->warpPaths;
115  this->continuousInputDataBuffer = ptr->continuousInputDataBuffer;
116  this->numTemplates = ptr->numTemplates;
117  this->rejectionMode = ptr->rejectionMode;
118  this->useSmoothing = ptr->useSmoothing;
119  this->useZNormalisation = ptr->useZNormalisation;
120  this->constrainZNorm = ptr->constrainZNorm;
121  this->constrainWarpingPath = ptr->constrainWarpingPath;
122  this->trimTrainingData = ptr->trimTrainingData;
123  this->zNormConstrainThreshold = ptr->zNormConstrainThreshold;
124  this->radius = ptr->radius;
125  this->offsetUsingFirstSample = ptr->offsetUsingFirstSample;
126  this->trimThreshold = ptr->trimThreshold;
127  this->maximumTrimPercentage = ptr->maximumTrimPercentage;
128  this->smoothingFactor = ptr->smoothingFactor;
129  this->distanceMethod = ptr->distanceMethod;
130  this->rejectionMode = ptr->rejectionMode;
131  this->averageTemplateLength = ptr->averageTemplateLength;
132 
133  //Copy the classifier variables
134  return copyBaseVariables( classifier );
135  }
136  return false;
137 }
138 
140 bool DTW::train_(TimeSeriesClassificationData &labelledTrainingData){
141 
142  UINT bestIndex = 0;
143 
144  //Cleanup Memory
145  templatesBuffer.clear();
146  classLabels.clear();
147  trained = false;
148  continuousInputDataBuffer.clear();
149 
150  if( trimTrainingData ){
151  TimeSeriesClassificationSampleTrimmer timeSeriesTrimmer(trimThreshold,maximumTrimPercentage);
153  tempData.setNumDimensions( labelledTrainingData.getNumDimensions() );
154 
155  for(UINT i=0; i<labelledTrainingData.getNumSamples(); i++){
156  if( timeSeriesTrimmer.trimTimeSeries( labelledTrainingData[i] ) ){
157  tempData.addSample(labelledTrainingData[i].getClassLabel(), labelledTrainingData[i].getData());
158  }else{
159  trainingLog << "Removing training sample " << i << " from the dataset as it could not be trimmed!" << endl;
160  }
161  }
162  //Overwrite the original training data with the trimmed dataset
163  labelledTrainingData = tempData;
164  }
165 
166  if( labelledTrainingData.getNumSamples() == 0 ){
167  errorLog << "train_(TimeSeriesClassificationData &labelledTrainingData) - Can't train model as there are no samples in training data!" << endl;
168  return false;
169  }
170 
171  //Assign
172  numClasses = labelledTrainingData.getNumClasses();
173  numTemplates = labelledTrainingData.getNumClasses();
174  numInputDimensions = labelledTrainingData.getNumDimensions();
175  templatesBuffer.resize( numClasses );
176  classLabels.resize( numClasses );
177  nullRejectionThresholds.resize( numClasses );
178  averageTemplateLength = 0;
179 
180  //Need to copy the labelled training data incase we need to scale it or znorm it
181  TimeSeriesClassificationData trainingData( labelledTrainingData );
182 
183  //Perform any scaling or normalisation
184  ranges = trainingData.getRanges();
185  if( useScaling ) scaleData( trainingData );
186  if( useZNormalisation ) znormData( trainingData );
187 
188  //For each class, run a one-to-one DTW and find the template the best describes the data
189  for(UINT k=0; k<numTemplates; k++){
190  //Get the class label for the cth class
191  UINT classLabel = trainingData.getClassTracker()[k].classLabel;
192  TimeSeriesClassificationData classData = trainingData.getClassData( classLabel );
193  UINT numExamples = classData.getNumSamples();
194  bestIndex = 0;
195 
196  //Set the class label of this template
197  templatesBuffer[k].classLabel = classLabel;
198 
199  //Set the kth class label
200  classLabels[k] = classLabel;
201 
202  trainingLog << "Training Template: " << k << " Class: " << classLabel << endl;
203 
204  //Check to make sure we actually have some training examples
205  if( numExamples < 1 ){
206  errorLog << "train_(TimeSeriesClassificationData &labelledTrainingData) - Can not train model: Num of Example is < 1! Class: " << classLabel << ". Turn off null rejection if you want to use DTW with only 1 training sample per class." << endl;
207  return false;
208  }
209 
210  if( numExamples == 1 && useNullRejection ){
211  errorLog << "train_(TimeSeriesClassificationData &labelledTrainingData) - Can not train model as there is only 1 example in class: " << classLabel << ". Turn off null rejection if you want to use DTW with only 1 training sample per class." << endl;
212  return false;
213  }
214 
215  if( numExamples == 1 ){//If we have just one training example then we have to use it as the template
216  bestIndex = 0;
217  nullRejectionThresholds[k] = 0.0;//TODO-We need a better way of calculating this!
218  }else{
219  //Search for the best training example for this class
220  if( !train_NDDTW(classData,templatesBuffer[k],bestIndex) ){
221  errorLog << "train_(LabelledTimeSeriesClassificationData &labelledTrainingData) - Failed to train template for class with label: " << classLabel << endl;
222  return false;
223  }
224  }
225 
226  //Add the template with the best index to the buffer
227  int trainingMethod = 0;
228  if(useSmoothing) trainingMethod = 1;
229 
230  switch (trainingMethod) {
231  case(0)://Standard Training
232  templatesBuffer[k].timeSeries = classData[bestIndex].getData();
233  break;
234  case(1)://Training using Smoothing
235  //Smooth the data, reducing its size by a factor set by smoothFactor
236  smoothData(classData[ bestIndex ].getData(),smoothingFactor,templatesBuffer[k].timeSeries);
237  break;
238  default:
239  cout<<"Can not train model: Unknown training method \n";
240  return false;
241  break;
242  }
243 
244  if( offsetUsingFirstSample ){
245  offsetTimeseries( templatesBuffer[k].timeSeries );
246  }
247 
248  //Add the average length of the training examples for this template to the overall averageTemplateLength
249  averageTemplateLength += templatesBuffer[k].averageTemplateLength;
250  }
251 
252  //Flag that the models have been trained
253  trained = true;
254  averageTemplateLength = averageTemplateLength/numTemplates;
255 
256  //Recompute the null rejection thresholds
257  recomputeNullRejectionThresholds();
258 
259  //Resize the prediction results to make sure it is setup for realtime prediction
260  continuousInputDataBuffer.clear();
261  continuousInputDataBuffer.resize(averageTemplateLength,vector<double>(numInputDimensions,0));
262  classLikelihoods.resize(numTemplates,DEFAULT_NULL_LIKELIHOOD_VALUE);
263  classDistances.resize(numTemplates,0);
264  predictedClassLabel = GRT_DEFAULT_NULL_CLASS_LABEL;
265  maxLikelihood = DEFAULT_NULL_LIKELIHOOD_VALUE;
266 
267  //Training complete
268  return true;
269 }
270 
271 bool DTW::train_NDDTW(TimeSeriesClassificationData &trainingData,DTWTemplate &dtwTemplate,UINT &bestIndex){
272 
273  UINT numExamples = trainingData.getNumSamples();
274  VectorDouble results(numExamples,0.0);
275  MatrixDouble distanceResults(numExamples,numExamples);
276  dtwTemplate.averageTemplateLength = 0;
277 
278  for(UINT m=0; m<numExamples; m++){
279 
280  MatrixDouble templateA; //The m'th template
281  MatrixDouble templateB; //The n'th template
282  dtwTemplate.averageTemplateLength += trainingData[m].getLength();
283 
284  //Smooth the data if required
285  if( useSmoothing ) smoothData(trainingData[m].getData(),smoothingFactor,templateA);
286  else templateA = trainingData[m].getData();
287 
288  if( offsetUsingFirstSample ){
289  offsetTimeseries(templateA);
290  }
291 
292  for(UINT n=0; n<numExamples; n++){
293  if(m!=n){
294  //Smooth the data if required
295  if( useSmoothing ) smoothData(trainingData[n].getData(),smoothingFactor,templateB);
296  else templateB = trainingData[n].getData();
297 
298  if( offsetUsingFirstSample ){
299  offsetTimeseries(templateB);
300  }
301 
302  //Compute the distance between the two time series
303  MatrixDouble distanceMatrix(templateA.getNumRows(),templateB.getNumRows());
304  vector< IndexDist > warpPath;
305  double dist = computeDistance(templateA,templateB,distanceMatrix,warpPath);
306 
307  trainingLog << "Template: " << m << " Timeseries: " << n << " Dist: " << dist << endl;
308 
309  //Update the results values
310  distanceResults[m][n] = dist;
311  results[m] += dist;
312  }else distanceResults[m][n] = 0; //The distance is zero because the two timeseries are the same
313  }
314  }
315 
316  for(UINT m=0; m<numExamples; m++) results[m]/=(numExamples-1);
317  //Find the best average result, this is the result with the minimum value
318  bestIndex = 0;
319  double bestAverage = results[0];
320  for(UINT m=1; m<numExamples; m++){
321  if( results[m] < bestAverage ){
322  bestAverage = results[m];
323  bestIndex = m;
324  }
325  }
326 
327  if( numExamples > 2 ){
328 
329  //Work out the threshold value for the best template
330  dtwTemplate.trainingMu = results[bestIndex];
331  dtwTemplate.trainingSigma = 0.0;
332 
333  for(UINT n=0; n<numExamples; n++){
334  if(n!=bestIndex){
335  dtwTemplate.trainingSigma += SQR( distanceResults[ bestIndex ][n] - dtwTemplate.trainingMu );
336  }
337  }
338  dtwTemplate.trainingSigma = sqrt( dtwTemplate.trainingSigma / double(numExamples-2) );
339  }else{
340  warningLog << "_train_NDDTW(TimeSeriesClassificationData &trainingData,DTWTemplate &dtwTemplate,UINT &bestIndex - There are not enough examples to compute the trainingMu and trainingSigma for the template for class " << dtwTemplate.classLabel << endl;
341  dtwTemplate.trainingMu = 0.0;
342  dtwTemplate.trainingSigma = 0.0;
343  }
344 
345  //Set the average length of the training examples
346  dtwTemplate.averageTemplateLength = (UINT) (dtwTemplate.averageTemplateLength/double(numExamples));
347 
348  trainingLog << "AverageTemplateLength: " << dtwTemplate.averageTemplateLength << endl;
349 
350  //Flag that the training was successfull
351  return true;
352 }
353 
354 
355 bool DTW::predict_(MatrixDouble &inputTimeSeries){
356 
357  if( !trained ){
358  errorLog << "predict_(MatrixDouble &inputTimeSeries) - The DTW templates have not been trained!" << endl;
359  return false;
360  }
361 
362  if( classLikelihoods.size() != numTemplates ) classLikelihoods.resize(numTemplates);
363  if( classDistances.size() != numTemplates ) classDistances.resize(numTemplates);
364 
365  predictedClassLabel = 0;
366  maxLikelihood = DEFAULT_NULL_LIKELIHOOD_VALUE;
367  for(UINT k=0; k<classLikelihoods.size(); k++){
368  classLikelihoods[k] = 0;
369  classDistances[k] = DEFAULT_NULL_LIKELIHOOD_VALUE;
370  }
371 
372  if( numInputDimensions != inputTimeSeries.getNumCols() ){
373  errorLog << "predict_(MatrixDouble &inputTimeSeries) - The number of features in the model (" << numInputDimensions << ") do not match that of the input time series (" << inputTimeSeries.getNumCols() << ")" << endl;
374  return false;
375  }
376 
377  //Perform any preprocessing if requried
378  MatrixDouble *timeSeriesPtr = &inputTimeSeries;
379  MatrixDouble processedTimeSeries;
380  MatrixDouble tempMatrix;
381  if(useScaling){
382  scaleData(*timeSeriesPtr,processedTimeSeries);
383  timeSeriesPtr = &processedTimeSeries;
384  }
385 
386  //Normalize the data if needed
387  if( useZNormalisation ){
388  znormData(*timeSeriesPtr,processedTimeSeries);
389  timeSeriesPtr = &processedTimeSeries;
390  }
391 
392  //Smooth the data if required
393  if( useSmoothing ){
394  smoothData(*timeSeriesPtr,smoothingFactor,tempMatrix);
395  timeSeriesPtr = &tempMatrix;
396  }
397 
398  //Offset the timeseries if required
399  if( offsetUsingFirstSample ){
400  offsetTimeseries( *timeSeriesPtr );
401  }
402 
403  //Make the prediction by finding the closest template
404  double sum = 0;
405  if( distanceMatrices.size() != numTemplates ) distanceMatrices.resize( numTemplates );
406  if( warpPaths.size() != numTemplates ) warpPaths.resize( numTemplates );
407 
408  //Test the timeSeries against all the templates in the timeSeries buffer
409  for(UINT k=0; k<numTemplates; k++){
410  //Perform DTW
411  classDistances[k] = computeDistance(templatesBuffer[k].timeSeries,*timeSeriesPtr,distanceMatrices[k],warpPaths[k]);
412  classLikelihoods[k] = 1.0 / classDistances[k];
413  sum += classLikelihoods[k];
414  }
415 
416  //See which gave the min distance
417  UINT closestTemplateIndex = 0;
418  bestDistance = classDistances[0];
419  for(UINT k=1; k<numTemplates; k++){
420  if( classDistances[k] < bestDistance ){
421  bestDistance = classDistances[k];
422  closestTemplateIndex = k;
423  }
424  }
425 
426  //Normalize the class likelihoods and check which class has the maximum likelihood
427  UINT maxLikelihoodIndex = 0;
428  maxLikelihood = 0;
429  if( sum > 0 ){
430  for(UINT k=0; k<numTemplates; k++){
431  classLikelihoods[k] /= sum;
432  if( classLikelihoods[k] > maxLikelihood ){
433  maxLikelihood = classLikelihoods[k];
434  maxLikelihoodIndex = k;
435  }
436  }
437  }
438 
439  if( useNullRejection ){
440 
441  switch( rejectionMode ){
442  case TEMPLATE_THRESHOLDS:
443  if( bestDistance <= nullRejectionThresholds[ closestTemplateIndex ] ) predictedClassLabel = templatesBuffer[ closestTemplateIndex ].classLabel;
444  else predictedClassLabel = GRT_DEFAULT_NULL_CLASS_LABEL;
445  break;
446  case CLASS_LIKELIHOODS:
447  if( maxLikelihood >= 0.99 ) predictedClassLabel = templatesBuffer[ maxLikelihoodIndex ].classLabel;
448  else predictedClassLabel = GRT_DEFAULT_NULL_CLASS_LABEL;
449  break;
450  case THRESHOLDS_AND_LIKELIHOODS:
451  if( bestDistance <= nullRejectionThresholds[ closestTemplateIndex ] && maxLikelihood >= 0.99 )
452  predictedClassLabel = templatesBuffer[ closestTemplateIndex ].classLabel;
453  else predictedClassLabel = GRT_DEFAULT_NULL_CLASS_LABEL;
454  break;
455  default:
456  errorLog << "predict_(MatrixDouble &timeSeries) - Unknown RejectionMode!" << endl;
457  return false;
458  break;
459  }
460 
461  }else predictedClassLabel = templatesBuffer[ closestTemplateIndex ].classLabel;
462 
463  return true;
464 }
465 
466 bool DTW::predict_(VectorDouble &inputVector){
467 
468  if( !trained ){
469  errorLog << "predict_(VectorDouble &inputVector) - The model has not been trained!" << endl;
470  return false;
471  }
472  predictedClassLabel = 0;
473  maxLikelihood = DEFAULT_NULL_LIKELIHOOD_VALUE;
474  std::fill(classLikelihoods.begin(),classLikelihoods.end(),DEFAULT_NULL_LIKELIHOOD_VALUE);
475  std::fill(classDistances.begin(),classDistances.end(),0);
476 
477  if( numInputDimensions != inputVector.size() ){
478  errorLog << "predict_(VectorDouble &inputVector) - The number of features in the model " << numInputDimensions << " does not match that of the input vector " << inputVector.size() << endl;
479  return false;
480  }
481 
482  //Add the new input to the circular buffer
483  continuousInputDataBuffer.push_back( inputVector );
484 
485  if( continuousInputDataBuffer.getNumValuesInBuffer() < averageTemplateLength ){
486  //We haven't got enough samples yet so can't do the prediction
487  return true;
488  }
489 
490  //Copy the data into a temporary matrix
491  const UINT M = continuousInputDataBuffer.getSize();
492  const UINT N = numInputDimensions;
493  MatrixDouble predictionTimeSeries(M,N);
494  for(UINT i=0; i<M; i++){
495  for(UINT j=0; j<N; j++){
496  predictionTimeSeries[i][j] = continuousInputDataBuffer[i][j];
497  }
498  }
499 
500  //Run the prediction
501  return predict( predictionTimeSeries );
502 
503 }
504 
505 bool DTW::reset(){
506  continuousInputDataBuffer.clear();
507  if( trained ){
508  continuousInputDataBuffer.resize(averageTemplateLength,vector<double>(numInputDimensions,0));
509  recomputeNullRejectionThresholds();
510  }
511  return true;
512 }
513 
514 bool DTW::clear(){
515 
516  //Clear the Classifier variables
517  Classifier::clear();
518 
519  //Clear the DTW model
520  templatesBuffer.clear();
521  distanceMatrices.clear();
522  warpPaths.clear();
523  continuousInputDataBuffer.clear();
524 
525  return true;
526 }
527 
528 bool DTW::recomputeNullRejectionThresholds(){
529 
530  if(!trained) return false;
531 
532  //Copy the null rejection thresholds into one buffer so they can easily be accessed from the base class
533  nullRejectionThresholds.resize(numTemplates);
534 
535  for(UINT k=0; k<numTemplates; k++){
536  //The threshold is set as the mean distance plus gamma standard deviations
537  nullRejectionThresholds[k] = templatesBuffer[k].trainingMu + (templatesBuffer[k].trainingSigma * nullRejectionCoeff);
538  }
539 
540  return true;
541 }
542 
543 bool DTW::setModels( vector< DTWTemplate > newTemplates ){
544 
545  if( newTemplates.size() == templatesBuffer.size() ){
546  templatesBuffer = newTemplates;
547  //Make sure the class labels have not changed
548  classLabels.resize( templatesBuffer.size() );
549  for(UINT i=0; i<templatesBuffer.size(); i++){
550  classLabels[i] = templatesBuffer[i].classLabel;
551  }
552  return true;
553  }
554  return false;
555 }
556 
558 
559 double DTW::computeDistance(MatrixDouble &timeSeriesA,MatrixDouble &timeSeriesB,MatrixDouble &distanceMatrix,vector< IndexDist > &warpPath){
560 
561  const int M = timeSeriesA.getNumRows();
562  const int N = timeSeriesB.getNumRows();
563  const int C = timeSeriesA.getNumCols();
564  int i,j,k,index = 0;
565  double totalDist,v,normFactor = 0.;
566 
567  warpPath.clear();
568  if( int(distanceMatrix.getNumRows()) != M || int(distanceMatrix.getNumCols()) != N ){
569  distanceMatrix.resize(M, N);
570  }
571 
572  switch (distanceMethod) {
573  case (ABSOLUTE_DIST):
574  for(i=0; i<M; i++){
575  for(j=0; j<N; j++){
576  distanceMatrix[i][j] = 0.0;
577  for(k=0; k< C; k++){
578  distanceMatrix[i][j] += fabs(timeSeriesA[i][k]-timeSeriesB[j][k]);
579  }
580  }
581  }
582  break;
583  case (EUCLIDEAN_DIST):
584  //Calculate Euclidean Distance for all possible values
585  for(i=0; i<M; i++){
586  for(j=0; j<N; j++){
587  distanceMatrix[i][j] = 0.0;
588  for(k=0; k< C; k++){
589  distanceMatrix[i][j] += SQR( timeSeriesA[i][k]-timeSeriesB[j][k] );
590  }
591  distanceMatrix[i][j] = sqrt( distanceMatrix[i][j] );
592  }
593  }
594  break;
595  case (NORM_ABSOLUTE_DIST):
596  for(i=0; i<M; i++){
597  for(j=0; j<N; j++){
598  distanceMatrix[i][j] = 0.0;
599  for(k=0; k< C; k++){
600  distanceMatrix[i][j] += fabs(timeSeriesA[i][k]-timeSeriesB[j][k]);
601  }
602  distanceMatrix[i][j]/=N;
603  }
604  }
605  break;
606  default:
607  errorLog<<"ERROR: Unknown distance method: "<<distanceMethod<<endl;
608  return -1;
609  break;
610  }
611 
612  //Run the recursive search function to build the cost matrix
613  double distance = sqrt( d(M-1,N-1,distanceMatrix,M,N) );
614 
615  if( grt_isinf(distance) || grt_isnan(distance) ){
616  warningLog << "DTW computeDistance(...) - Distance Matrix Values are INF!" << endl;
617  return INFINITY;
618  }
619 
620  //cout << "DIST: " << distance << endl;
621 
622  //The distMatrix values are negative so make them positive
623  for(i=0; i<M; i++){
624  for(j=0; j<N; j++){
625  distanceMatrix[i][j] = fabs( distanceMatrix[i][j] );
626  }
627  }
628 
629  //Now Create the Warp Path through the cost matrix, starting at the end
630  i=M-1;
631  j=N-1;
632  totalDist = distanceMatrix[i][j];
633  warpPath.push_back( IndexDist(i,j,distanceMatrix[i][j]) );
634 
635  //Use dynamic programming to navigate through the cost matrix until [0][0] has been reached
636  normFactor = 1;
637  while( true ) {
638  if( i==0 && j==0 ) break;
639  if( i==0 ){ j--; }
640  else{
641  if( j==0 ) i--;
642  else{
643  //Find the minimum cell to move to
644  v = numeric_limits<double>::max();
645  index = 0;
646  if( distanceMatrix[i-1][j] < v ){ v = distanceMatrix[i-1][j]; index = 1; }
647  if( distanceMatrix[i][j-1] < v ){ v = distanceMatrix[i][j-1]; index = 2; }
648  if( distanceMatrix[i-1][j-1] <= v ){ index = 3; }
649  switch(index){
650  case(1):
651  i--;
652  break;
653  case(2):
654  j--;
655  break;
656  case(3):
657  i--;
658  j--;
659  break;
660  default:
661  warningLog << "DTW computeDistance(...) - Could not compute a warping path for the input matrix! Dist: " << distanceMatrix[i-1][j] << " i: " << i << " j: " << j << endl;
662  return INFINITY;
663  break;
664  }
665  }
666  }
667  normFactor++;
668  totalDist += distanceMatrix[i][j];
669  warpPath.push_back( IndexDist(i,j,distanceMatrix[i][j]) );
670  }
671 
672  return totalDist/normFactor;
673 }
674 
675 double DTW::d(int m,int n,MatrixDouble &distanceMatrix,const int M,const int N){
676  double dist = 0;
677  //The following is based on Matlab code by Eamonn Keogh and Michael Pazzani
678 
679  //If this cell is NAN then it has already been flagged as unreachable
680  if( grt_isnan( distanceMatrix[m][n] ) ){
681  return NAN;
682  }
683 
684  if( constrainWarpingPath ){
685  double r = ceil( min(M,N)*radius );
686  //Test to see if the current cell is outside of the warping window
687  if( fabs( n-((N-1)/((M-1)/double(m))) ) > r ){
688  if( n-((N-1)/((M-1)/double(m))) > 0 ){
689  for(int i=0; i<m; i++){
690  for(int j=n; j<N; j++){
691  distanceMatrix[i][j] = NAN;
692  }
693  }
694  }else{
695  for(int i=m; i<M; i++){
696  for(int j=0; j<n; j++){
697  distanceMatrix[i][j] = NAN;
698  }
699  }
700  }
701  return NAN;
702  }
703  }
704 
705  //If this cell contains a negative value then it has already been searched
706  //The cost is therefore the absolute value of the negative value so return it
707  if( distanceMatrix[m][n] < 0 ){
708  dist = fabs( distanceMatrix[m][n] );
709  return dist;
710  }
711 
712  //Case 1: A warping path has reached the end
713  //Return the contribution of distance
714  //Negate the value, to record the fact that this cell has been visited
715  //End of recursion
716  if( m == 0 && n == 0 ){
717  dist = distanceMatrix[0][0];
718  distanceMatrix[0][0] = -distanceMatrix[0][0];
719  return dist;
720  }
721 
722  //Case 2: we are somewhere in the top row of the matrix
723  //Only need to consider moving left
724  if( m == 0 ){
725  double contribDist = d(m,n-1,distanceMatrix,M,N);
726 
727  dist = distanceMatrix[m][n] + contribDist;
728 
729  distanceMatrix[m][n] = -dist;
730  return dist;
731  }else{
732  //Case 3: we are somewhere in the left column of the matrix
733  //Only need to consider moving down
734  if ( n == 0) {
735  double contribDist = d(m-1,n,distanceMatrix,M,N);
736 
737  dist = distanceMatrix[m][n] + contribDist;
738 
739  distanceMatrix[m][n] = -dist;
740  return dist;
741  }else{
742  //Case 4: We are somewhere away from the edges so consider moving in the three main directions
743  double contribDist1 = d(m-1,n-1,distanceMatrix,M,N);
744  double contribDist2 = d(m-1,n,distanceMatrix,M,N);
745  double contribDist3 = d(m,n-1,distanceMatrix,M,N);
746  double minValue = numeric_limits<double>::max();
747  int index = 0;
748  if( contribDist1 < minValue ){ minValue = contribDist1; index = 1; }
749  if( contribDist2 < minValue ){ minValue = contribDist2; index = 2; }
750  if( contribDist3 < minValue ){ minValue = contribDist3; index = 3; }
751 
752  switch ( index ) {
753  case 1:
754  dist = distanceMatrix[m][n] + minValue;
755  break;
756  case 2:
757  dist = distanceMatrix[m][n] + minValue;
758  break;
759  case 3:
760  dist = distanceMatrix[m][n] + minValue;
761  break;
762 
763  default:
764  break;
765  }
766 
767  distanceMatrix[m][n] = -dist; //Negate the value to record that it has been visited
768  return dist;
769  }
770  }
771 
772  //This should not happen!
773  return dist;
774 }
775 
776 inline double DTW::MIN_(double a,double b, double c){
777  double v = a;
778  if(b<v) v = b;
779  if(c<v) v = c;
780  return v;
781 }
782 
783 
785 
786 void DTW::scaleData(TimeSeriesClassificationData &trainingData){
787 
788  //Scale the data using the min and max values
789  for(UINT i=0; i<trainingData.getNumSamples(); i++){
790  scaleData( trainingData[i].getData(), trainingData[i].getData() );
791  }
792 
793 }
794 
795 void DTW::scaleData(MatrixDouble &data,MatrixDouble &scaledData){
796 
797  const UINT R = data.getNumRows();
798  const UINT C = data.getNumCols();
799 
800  if( scaledData.getNumRows() != R || scaledData.getNumCols() != C ){
801  scaledData.resize(R, C);
802  }
803 
804  //Scale the data using the min and max values
805  for(UINT i=0; i<R; i++)
806  for(UINT j=0; j<C; j++)
807  scaledData[i][j] = scale(data[i][j],ranges[j].minValue,ranges[j].maxValue,0.0,1.0);
808 
809 }
810 
811 void DTW::znormData(TimeSeriesClassificationData &trainingData){
812 
813  for(UINT i=0; i<trainingData.getNumSamples(); i++){
814  znormData( trainingData[i].getData(), trainingData[i].getData() );
815  }
816 
817 }
818 
819 void DTW::znormData(MatrixDouble &data,MatrixDouble &normData){
820 
821  const UINT R = data.getNumRows();
822  const UINT C = data.getNumCols();
823 
824  if( normData.getNumRows() != R || normData.getNumCols() != C ){
825  normData.resize(R,C);
826  }
827 
828  for(UINT j=0; j<C; j++){
829  double mean = 0.0;
830  double stdDev = 0.0;
831 
832  //Calculate Mean
833  for(UINT i=0; i<R; i++) mean += data[i][j];
834  mean /= double(R);
835 
836  //Calculate Std Dev
837  for(UINT i=0; i<R; i++)
838  stdDev += SQR(data[i][j]-mean);
839  stdDev = sqrt( stdDev / (R - 1.0) );
840 
841  if(constrainZNorm && stdDev < 0.01){
842  //Normalize the data to 0 mean
843  for(UINT i=0; i<R; i++)
844  normData[i][j] = (data[i][j] - mean);
845  }else{
846  //Normalize the data to 0 mean and standard deviation of 1
847  for(UINT i=0; i<R; i++)
848  normData[i][j] = (data[i][j] - mean) / stdDev;
849  }
850  }
851 }
852 
853 void DTW::smoothData(VectorDouble &data,UINT smoothFactor,VectorDouble &resultsData){
854 
855  const UINT M = (UINT)data.size();
856  const UINT N = (UINT) floor(double(M)/double(smoothFactor));
857  resultsData.resize(N,0);
858  for(UINT i=0; i<N; i++) resultsData[i]=0.0;
859 
860  if(smoothFactor==1 || M<smoothFactor){
861  resultsData = data;
862  return;
863  }
864 
865  for(UINT i=0; i<N; i++){
866  double mean = 0.0;
867  UINT index = i*smoothFactor;
868  for(UINT x=0; x<smoothFactor; x++){
869  mean += data[index+x];
870  }
871  resultsData[i] = mean/smoothFactor;
872  }
873  //Add on the data that does not fit into the window
874  if(M%smoothFactor!=0.0){
875  double mean = 0.0;
876  for(UINT i=N*smoothFactor; i<M; i++) mean += data[i];
877  mean/=M-(N*smoothFactor);
878  //Add one to the end of the vector
879  VectorDouble tempVector(N+1);
880  for(UINT i=0; i<N; i++) tempVector[i] = resultsData[i];
881  tempVector[N] = mean;
882  resultsData = tempVector;
883  }
884 
885 }
886 
887 void DTW::smoothData(MatrixDouble &data,UINT smoothFactor,MatrixDouble &resultsData){
888 
889  const UINT M = data.getNumRows();
890  const UINT C = data.getNumCols();
891  const UINT N = (UINT) floor(double(M)/double(smoothFactor));
892  resultsData.resize(N,C);
893 
894  if(smoothFactor==1 || M<smoothFactor){
895  resultsData = data;
896  return;
897  }
898 
899  for(UINT i=0; i<N; i++){
900  for(UINT j=0; j<C; j++){
901  double mean = 0.0;
902  int index = i*smoothFactor;
903  for(UINT x=0; x<smoothFactor; x++){
904  mean += data[index+x][j];
905  }
906  resultsData[i][j] = mean/smoothFactor;
907  }
908  }
909 
910  //Add on the data that does not fit into the window
911  if(M%smoothFactor!=0.0){
912  VectorDouble mean(C,0.0);
913  for(UINT j=0; j<C; j++){
914  for(UINT i=N*smoothFactor; i<M; i++) mean[j] += data[i][j];
915  mean[j]/=M-(N*smoothFactor);
916  }
917 
918  //Add one row to the end of the Matrix
919  MatrixDouble tempMatrix(N+1,C);
920 
921  for(UINT i=0; i<N; i++)
922  for(UINT j=0; j<C; j++)
923  tempMatrix[i][j] = resultsData[i][j];
924 
925  for(UINT j=0; j<C; j++) tempMatrix[N][j] = mean[j];
926  resultsData = tempMatrix;
927  }
928 
929 }
930 
932 
933 bool DTW::saveModelToFile( fstream &file ) const{
934 
935  if(!file.is_open()){
936  errorLog << "saveModelToFile( string fileName ) - Could not open file to save data" << endl;
937  return false;
938  }
939 
940  file << "GRT_DTW_Model_File_V2.0" <<endl;
941 
942  //Write the classifier settings to the file
943  if( !Classifier::saveBaseSettingsToFile(file) ){
944  errorLog <<"saveModelToFile(fstream &file) - Failed to save classifier base settings to file!" << endl;
945  return false;
946  }
947 
948  file << "DistanceMethod: ";
949  switch(distanceMethod){
950  case(ABSOLUTE_DIST):
951  file <<ABSOLUTE_DIST<<endl;
952  break;
953  case(EUCLIDEAN_DIST):
954  file <<EUCLIDEAN_DIST<<endl;
955  break;
956  default:
957  file <<ABSOLUTE_DIST<<endl;
958  break;
959  }
960  file << "UseSmoothing: "<<useSmoothing<<endl;
961  file << "SmoothingFactor: "<<smoothingFactor<<endl;
962  file << "UseZNormalisation: "<<useZNormalisation<<endl;
963  file << "OffsetUsingFirstSample: " << offsetUsingFirstSample << endl;
964  file << "ConstrainWarpingPath: " << constrainWarpingPath << endl;
965  file << "Radius: " << radius << endl;
966  file << "RejectionMode: " << rejectionMode<< endl;
967 
968  if( trained ){
969  file << "NumberOfTemplates: " << numTemplates << endl;
970  file << "OverallAverageTemplateLength: " << averageTemplateLength << endl;
971  //Save each template
972  for(UINT i=0; i<numTemplates; i++){
973  file << "***************TEMPLATE***************" << endl;
974  file << "Template: " << i+1 << endl;
975  file << "ClassLabel: " << templatesBuffer[i].classLabel << endl;
976  file << "TimeSeriesLength: " << templatesBuffer[i].timeSeries.getNumRows() << endl;
977  file << "TemplateThreshold: " << nullRejectionThresholds[i] << endl;
978  file << "TrainingMu: " << templatesBuffer[i].trainingMu << endl;
979  file << "TrainingSigma: " << templatesBuffer[i].trainingSigma << endl;
980  file << "AverageTemplateLength: " << templatesBuffer[i].averageTemplateLength << endl;
981  file << "TimeSeries: " << endl;
982  for(UINT k=0; k<templatesBuffer[i].timeSeries.getNumRows(); k++){
983  for(UINT j=0; j<templatesBuffer[i].timeSeries.getNumCols(); j++){
984  file << templatesBuffer[i].timeSeries[k][j] << "\t";
985  }
986  file << endl;
987  }
988  }
989  }
990 
991  return true;
992 }
993 
994 bool DTW::loadModelFromFile( fstream &file ){
995 
996  std::string word;
997  UINT timeSeriesLength;
998  UINT ts;
999 
1000  if(!file.is_open())
1001  {
1002  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to open file!" << endl;
1003  return false;
1004  }
1005 
1006  file >> word;
1007 
1008  //Check to see if we should load a legacy file
1009  if( word == "GRT_DTW_Model_File_V1.0" ){
1010  return loadLegacyModelFromFile( file );
1011  }
1012 
1013  //Check to make sure this is a file with the DTW File Format
1014  if(word != "GRT_DTW_Model_File_V2.0"){
1015  errorLog << "loadDTWModelFromFile( string fileName ) - Unknown file header!" << endl;
1016  return false;
1017  }
1018 
1019  //Load the base settings from the file
1020  if( !Classifier::loadBaseSettingsFromFile(file) ){
1021  errorLog << "loadModelFromFile(string filename) - Failed to load base settings from file!" << endl;
1022  return false;
1023  }
1024 
1025  //Check and load the Distance Method
1026  file >> word;
1027  if(word != "DistanceMethod:"){
1028  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find DistanceMethod!" << endl;
1029  return false;
1030  }
1031  file >> distanceMethod;
1032 
1033  //Check and load if Smoothing is used
1034  file >> word;
1035  if(word != "UseSmoothing:"){
1036  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find UseSmoothing!" << endl;
1037  return false;
1038  }
1039  file >> useSmoothing;
1040 
1041  //Check and load what the smoothing factor is
1042  file >> word;
1043  if(word != "SmoothingFactor:"){
1044  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find SmoothingFactor!" << endl;
1045  return false;
1046  }
1047  file >> smoothingFactor;
1048 
1049  //Check and load if ZNormalization is used
1050  file >> word;
1051  if(word != "UseZNormalisation:"){
1052  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find UseZNormalisation!" << endl;
1053  return false;
1054  }
1055  file >> useZNormalisation;
1056 
1057  //Check and load if OffsetUsingFirstSample is used
1058  file >> word;
1059  if(word != "OffsetUsingFirstSample:"){
1060  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find OffsetUsingFirstSample!" << endl;
1061  return false;
1062  }
1063  file >> offsetUsingFirstSample;
1064 
1065  //Check and load if ConstrainWarpingPath is used
1066  file >> word;
1067  if(word != "ConstrainWarpingPath:"){
1068  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find ConstrainWarpingPath!" << endl;
1069  return false;
1070  }
1071  file >> constrainWarpingPath;
1072 
1073  //Check and load if ZNormalization is used
1074  file >> word;
1075  if(word != "Radius:"){
1076  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find Radius!" << endl;
1077  return false;
1078  }
1079  file >> radius;
1080 
1081  //Check and load if Scaling is used
1082  file >> word;
1083  if(word != "RejectionMode:"){
1084  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find RejectionMode!" << endl;
1085  return false;
1086  }
1087  file >> rejectionMode;
1088 
1089  if( trained ){
1090 
1091  //Check and load the Number of Templates
1092  file >> word;
1093  if(word != "NumberOfTemplates:"){
1094  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find NumberOfTemplates!" << endl;
1095  return false;
1096  }
1097  file >> numTemplates;
1098 
1099  //Check and load the overall average template length
1100  file >> word;
1101  if(word != "OverallAverageTemplateLength:"){
1102  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find OverallAverageTemplateLength!" << endl;
1103  return false;
1104  }
1105  file >> averageTemplateLength;
1106 
1107  //Clean and reset the memory
1108  templatesBuffer.resize(numTemplates);
1109  classLabels.resize(numTemplates);
1110  nullRejectionThresholds.resize(numTemplates);
1111 
1112  //Load each template
1113  for(UINT i=0; i<numTemplates; i++){
1114  //Check we have the correct template
1115  file >> word;
1116  if( word != "***************TEMPLATE***************" ){
1117  clear();
1118  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find template header!" << endl;
1119  return false;
1120  }
1121 
1122  //Load the template number
1123  file >> word;
1124  if(word != "Template:"){
1125  clear();
1126  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find Template Number!" << endl;
1127  return false;
1128  }
1129 
1130  //Check the template number
1131  file >> ts;
1132  if(ts!=i+1){
1133  clear();
1134  errorLog << "loadDTWModelFromFile( string fileName ) - Invalid Template Number: " << ts << endl;
1135  return false;
1136  }
1137 
1138  //Get the class label of this template
1139  file >> word;
1140  if(word != "ClassLabel:"){
1141  clear();
1142  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find ClassLabel!" << endl;
1143  return false;
1144  }
1145  file >> templatesBuffer[i].classLabel;
1146  classLabels[i] = templatesBuffer[i].classLabel;
1147 
1148  //Get the time series length
1149  file >> word;
1150  if(word != "TimeSeriesLength:"){
1151  clear();
1152  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find TimeSeriesLength!" << endl;
1153  return false;
1154  }
1155  file >> timeSeriesLength;
1156 
1157  //Resize the buffers
1158  templatesBuffer[i].timeSeries.resize(timeSeriesLength,numInputDimensions);
1159 
1160  //Get the template threshold
1161  file >> word;
1162  if(word != "TemplateThreshold:"){
1163  clear();
1164  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find TemplateThreshold!" << endl;
1165  return false;
1166  }
1167  file >> nullRejectionThresholds[i];
1168 
1169  //Get the mu values
1170  file >> word;
1171  if(word != "TrainingMu:"){
1172  clear();
1173  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find TrainingMu!" << endl;
1174  return false;
1175  }
1176  file >> templatesBuffer[i].trainingMu;
1177 
1178  //Get the sigma values
1179  file >> word;
1180  if(word != "TrainingSigma:"){
1181  clear();
1182  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find TrainingSigma!" << endl;
1183  return false;
1184  }
1185  file >> templatesBuffer[i].trainingSigma;
1186 
1187  //Get the AverageTemplateLength value
1188  file >> word;
1189  if(word != "AverageTemplateLength:"){
1190  clear();
1191  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find AverageTemplateLength!" << endl;
1192  return false;
1193  }
1194  file >> templatesBuffer[i].averageTemplateLength;
1195 
1196  //Get the data
1197  file >> word;
1198  if(word != "TimeSeries:"){
1199  clear();
1200  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find template timeseries!" << endl;
1201  return false;
1202  }
1203  for(UINT k=0; k<timeSeriesLength; k++)
1204  for(UINT j=0; j<numInputDimensions; j++)
1205  file >> templatesBuffer[i].timeSeries[k][j];
1206  }
1207 
1208  //Resize the prediction results to make sure it is setup for realtime prediction
1209  continuousInputDataBuffer.clear();
1210  continuousInputDataBuffer.resize(averageTemplateLength,vector<double>(numInputDimensions,0));
1211  maxLikelihood = DEFAULT_NULL_LIKELIHOOD_VALUE;
1212  bestDistance = DEFAULT_NULL_DISTANCE_VALUE;
1213  classLikelihoods.resize(numClasses,DEFAULT_NULL_LIKELIHOOD_VALUE);
1214  classDistances.resize(numClasses,DEFAULT_NULL_DISTANCE_VALUE);
1215  }
1216 
1217  return true;
1218 }
1219 bool DTW::setRejectionMode(UINT rejectionMode){
1220  if( rejectionMode == TEMPLATE_THRESHOLDS || rejectionMode == CLASS_LIKELIHOODS || rejectionMode == THRESHOLDS_AND_LIKELIHOODS ){
1221  this->rejectionMode = rejectionMode;
1222  return true;
1223  }
1224  return false;
1225 }
1226 
1227 bool DTW::setOffsetTimeseriesUsingFirstSample(bool offsetUsingFirstSample){
1228  this->offsetUsingFirstSample = offsetUsingFirstSample;
1229  return true;
1230 }
1231 
1232 bool DTW::setContrainWarpingPath(bool constrain){
1233  this->constrainWarpingPath = constrain;
1234  return true;
1235 }
1236 
1237 bool DTW::setWarpingRadius(double radius){
1238  this->radius = radius;
1239  return true;
1240 }
1241 
1242 bool DTW::enableZNormalization(bool useZNormalisation,bool constrainZNorm){
1243  this->useZNormalisation = useZNormalisation;
1244  this->constrainZNorm = constrainZNorm;
1245  return true;
1246 }
1247 
1248 bool DTW::enableTrimTrainingData(bool trimTrainingData,double trimThreshold,double maximumTrimPercentage){
1249 
1250  if( trimThreshold < 0 || trimThreshold > 1 ){
1251  warningLog << "Failed to set trimTrainingData. The trimThreshold must be in the range of [0 1]" << endl;
1252  return false;
1253  }
1254  if( maximumTrimPercentage < 0 || maximumTrimPercentage > 100 ){
1255  warningLog << "Failed to set trimTrainingData. The maximumTrimPercentage must be a valid percentage in the range of [0 100]" << endl;
1256  return false;
1257  }
1258 
1259  this->trimTrainingData = trimTrainingData;
1260  this->trimThreshold = trimThreshold;
1261  this->maximumTrimPercentage = maximumTrimPercentage;
1262  return true;
1263 }
1264 
1265 void DTW::offsetTimeseries(MatrixDouble &timeseries){
1266  VectorDouble firstRow = timeseries.getRowVector(0);
1267  for(UINT i=0; i<timeseries.getNumRows(); i++){
1268  for(UINT j=0; j<timeseries.getNumCols(); j++){
1269  timeseries[i][j] -= firstRow[j];
1270  }
1271  }
1272 }
1273 
1274 bool DTW::loadLegacyModelFromFile( fstream &file ){
1275 
1276  string word;
1277  UINT timeSeriesLength;
1278  UINT ts;
1279 
1280  //Check and load the Number of Dimensions
1281  file >> word;
1282  if(word != "NumberOfDimensions:"){
1283  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find NumberOfDimensions!" << endl;
1284  return false;
1285  }
1286  file >> numInputDimensions;
1287 
1288  //Check and load the Number of Classes
1289  file >> word;
1290  if(word != "NumberOfClasses:"){
1291  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find NumberOfClasses!" << endl;
1292  return false;
1293  }
1294  file >> numClasses;
1295 
1296  //Check and load the Number of Templates
1297  file >> word;
1298  if(word != "NumberOfTemplates:"){
1299  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find NumberOfTemplates!" << endl;
1300  return false;
1301  }
1302  file >> numTemplates;
1303 
1304  //Check and load the Distance Method
1305  file >> word;
1306  if(word != "DistanceMethod:"){
1307  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find DistanceMethod!" << endl;
1308  return false;
1309  }
1310  file >> distanceMethod;
1311 
1312  //Check and load if UseNullRejection is used
1313  file >> word;
1314  if(word != "UseNullRejection:"){
1315  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find UseNullRejection!" << endl;
1316  return false;
1317  }
1318  file >> useNullRejection;
1319 
1320  //Check and load if Smoothing is used
1321  file >> word;
1322  if(word != "UseSmoothing:"){
1323  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find UseSmoothing!" << endl;
1324  return false;
1325  }
1326  file >> useSmoothing;
1327 
1328  //Check and load what the smoothing factor is
1329  file >> word;
1330  if(word != "SmoothingFactor:"){
1331  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find SmoothingFactor!" << endl;
1332  return false;
1333  }
1334  file >> smoothingFactor;
1335 
1336  //Check and load if Scaling is used
1337  file >> word;
1338  if(word != "UseScaling:"){
1339  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find UseScaling!" << endl;
1340  return false;
1341  }
1342  file >> useScaling;
1343 
1344  //Check and load if ZNormalization is used
1345  file >> word;
1346  if(word != "UseZNormalisation:"){
1347  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find UseZNormalisation!" << endl;
1348  return false;
1349  }
1350  file >> useZNormalisation;
1351 
1352  //Check and load if OffsetUsingFirstSample is used
1353  file >> word;
1354  if(word != "OffsetUsingFirstSample:"){
1355  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find OffsetUsingFirstSample!" << endl;
1356  return false;
1357  }
1358  file >> offsetUsingFirstSample;
1359 
1360  //Check and load if ConstrainWarpingPath is used
1361  file >> word;
1362  if(word != "ConstrainWarpingPath:"){
1363  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find ConstrainWarpingPath!" << endl;
1364  return false;
1365  }
1366  file >> constrainWarpingPath;
1367 
1368  //Check and load if ZNormalization is used
1369  file >> word;
1370  if(word != "Radius:"){
1371  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find Radius!" << endl;
1372  return false;
1373  }
1374  file >> radius;
1375 
1376  //Check and load if Scaling is used
1377  file >> word;
1378  if(word != "RejectionMode:"){
1379  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find RejectionMode!" << endl;
1380  return false;
1381  }
1382  file >> rejectionMode;
1383 
1384  //Check and load gamma
1385  file >> word;
1386  if(word != "NullRejectionCoeff:"){
1387  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find NullRejectionCoeff!" << endl;
1388  return false;
1389  }
1390  file >> nullRejectionCoeff;
1391 
1392  //Check and load the overall average template length
1393  file >> word;
1394  if(word != "OverallAverageTemplateLength:"){
1395  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find OverallAverageTemplateLength!" << endl;
1396  return false;
1397  }
1398  file >> averageTemplateLength;
1399 
1400  //Clean and reset the memory
1401  templatesBuffer.resize(numTemplates);
1402  classLabels.resize(numTemplates);
1403  nullRejectionThresholds.resize(numTemplates);
1404 
1405  //Load each template
1406  for(UINT i=0; i<numTemplates; i++){
1407  //Check we have the correct template
1408  file >> word;
1409  while(word != "Template:"){
1410  file >> word;
1411  }
1412  file >> ts;
1413 
1414  //Check the template number
1415  if(ts!=i+1){
1416  numTemplates=0;
1417  trained = false;
1418  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find Invalid Template Number!" << endl;
1419  return false;
1420  }
1421 
1422  //Get the class label of this template
1423  file >> word;
1424  if(word != "ClassLabel:"){
1425  numTemplates=0;
1426  trained = false;
1427  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find ClassLabel!" << endl;
1428  return false;
1429  }
1430  file >> templatesBuffer[i].classLabel;
1431  classLabels[i] = templatesBuffer[i].classLabel;
1432 
1433  //Get the time series length
1434  file >> word;
1435  if(word != "TimeSeriesLength:"){
1436  numTemplates=0;
1437  trained = false;
1438  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find TimeSeriesLength!" << endl;
1439  return false;
1440  }
1441  file >> timeSeriesLength;
1442 
1443  //Resize the buffers
1444  templatesBuffer[i].timeSeries.resize(timeSeriesLength,numInputDimensions);
1445 
1446  //Get the template threshold
1447  file >> word;
1448  if(word != "TemplateThreshold:"){
1449  numTemplates=0;
1450  trained = false;
1451  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find TemplateThreshold!" << endl;
1452  return false;
1453  }
1454  file >> nullRejectionThresholds[i];
1455 
1456  //Get the mu values
1457  file >> word;
1458  if(word != "TrainingMu:"){
1459  numTemplates=0;
1460  trained = false;
1461  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find TrainingMu!" << endl;
1462  return false;
1463  }
1464  file >> templatesBuffer[i].trainingMu;
1465 
1466  //Get the sigma values
1467  file >> word;
1468  if(word != "TrainingSigma:"){
1469  numTemplates=0;
1470  trained = false;
1471  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find TrainingSigma!" << endl;
1472  return false;
1473  }
1474  file >> templatesBuffer[i].trainingSigma;
1475 
1476  //Get the AverageTemplateLength value
1477  file >> word;
1478  if(word != "AverageTemplateLength:"){
1479  numTemplates=0;
1480  trained = false;
1481  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find AverageTemplateLength!" << endl;
1482  return false;
1483  }
1484  file >> templatesBuffer[i].averageTemplateLength;
1485 
1486  //Get the data
1487  file >> word;
1488  if(word != "TimeSeries:"){
1489  numTemplates=0;
1490  trained = false;
1491  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find template timeseries!" << endl;
1492  return false;
1493  }
1494  for(UINT k=0; k<timeSeriesLength; k++)
1495  for(UINT j=0; j<numInputDimensions; j++)
1496  file >> templatesBuffer[i].timeSeries[k][j];
1497 
1498  //Check for the footer
1499  file >> word;
1500  if(word != "***************************"){
1501  numTemplates=0;
1502  numClasses = 0;
1503  numInputDimensions=0;
1504  trained = false;
1505  errorLog << "loadDTWModelFromFile( string fileName ) - Failed to find template footer!" << endl;
1506  return false;
1507  }
1508  }
1509 
1510  //Resize the prediction results to make sure it is setup for realtime prediction
1511  continuousInputDataBuffer.clear();
1512  continuousInputDataBuffer.resize(averageTemplateLength,vector<double>(numInputDimensions,0));
1513  maxLikelihood = DEFAULT_NULL_LIKELIHOOD_VALUE;
1514  bestDistance = DEFAULT_NULL_DISTANCE_VALUE;
1515  classLikelihoods.resize(numClasses,DEFAULT_NULL_LIKELIHOOD_VALUE);
1516  classDistances.resize(numClasses,DEFAULT_NULL_DISTANCE_VALUE);
1517 
1518  trained = true;
1519 
1520  return true;
1521 }
1522 
1523 } //End of namespace GRT
bool trimTimeSeries(TimeSeriesClassificationSample &timeSeries)
Definition: AdaBoost.cpp:25
bool push_back(const std::vector< T > &sample)
Definition: Matrix.h:390
unsigned int getNumCols() const
Definition: Matrix.h:538
bool setNumDimensions(const UINT numDimensions)
bool addSample(const UINT classLabel, const MatrixDouble &trainingSample)
Definition: DTW.h:91
This class implements Dynamic Time Warping. Dynamic Time Warping (DTW) is a powerful classifier that ...
vector< ClassTracker > getClassTracker() const
TimeSeriesClassificationData getClassData(const UINT classLabel) const
unsigned int getNumRows() const
Definition: Matrix.h:531
string getClassifierType() const
Definition: Classifier.cpp:159
std::vector< T > getRowVector(const unsigned int r) const
Definition: Matrix.h:173
virtual bool resize(const unsigned int r, const unsigned int c)
Definition: Matrix.h:234