Orfeo Toolbox  4.0
otbMeanShiftKernels.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: ORFEO Toolbox
4  Language: C++
5  Date: $Date$
6  Version: $Revision$
7 
8 
9  Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
10  See OTBCopyright.txt for details.
11 
12 
13  This software is distributed WITHOUT ANY WARRANTY; without even
14  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
15  PURPOSE. See the above copyright notices for more information.
16 
17 =========================================================================*/
18 #ifndef __otbMeanShiftKernels_h
19 #define __otbMeanShiftKernels_h
20 
21 #include <vector>
22 #include <algorithm> // for std::find
23 #include <cstring> // for strpbrk
24 
25 #include "itkNumericTraits.h"
26 #include "otbMath.h"
27 
28 // Existing kernels
33 
34 
35 /*Balls itk
36  *
37  *
38  */
39 
40 
41 /* SVM kernels */
42 /*
43  //
44  * */
45 
46 
47 namespace otb
48 {
53 {
54 public:
57 
59  {
60  this->SetName("custom");
61  }
62  virtual ~CustomKernelFunctor() {}
63 
64  // Deep copy operator
66  {
67  return new Self(*this);
68  }
69 
70  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& /*param*/) const
71  {
72  return (dot(x, x) - 2.0 * dot(x, y) + dot(y, y));
73  }
74 
75 protected:
77  : Superclass(copy)
78  {
79  *this = copy;
80  }
81 
83  {
85  return *this;
86  }
87 };
88 
93 {
94 public:
97 
99  {
100  this->SetName("invMultiQuadric");
101  this->SetValue<double>("const_coef", 1);
102  }
104 
105  // Deep copy operator
107  {
108  return new Self(*this);
109  }
110 
111  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
112  {
113  double mq = this->GetValue<double>("const_coef") + m_Custom(x, y, param);
114  if (mq == 0.0)
115  {
117  }
118  return 1.0 / sqrt(mq);
119  }
120 
121 protected:
123  : Superclass(copy)
124  {
125  *this = copy;
126  }
127 
129  {
131  return *this;
132  }
133 
134 private:
136 };
137 
142 {
143 public:
146 
148  {
149  this->SetName("KMod");
150  this->SetValue<double>("const_coef", 1);
151  }
152  virtual ~KModKernelFunctor() {}
153 
154  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
155  {
156  double mq = this->GetValue<double>("const_coef") + m_Custom(x, y, param);
157 
158  if (mq == 0.0)
159  {
161  }
162  return exp(param.gamma / mq) - 1.0;
163  }
164 
165  // Deep copy operator
167  {
168  return new Self(*this);
169  }
170 
171 protected:
172  KModKernelFunctor(const Self& copy)
173  : Superclass(copy)
174  {
175  *this = copy;
176  }
177 
179  {
181  return *this;
182  }
183 
184 private:
186 };
187 
192 {
193 public:
196 
198  {
199  this->SetName("SAM");
200  }
201  virtual ~SAMKernelFunctor() {}
202 
203  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& /*param*/) const
204  {
205  double den = dot(x, x) * dot(y, y);
206  if (den <= 0.)
207  {
208  return 0.0;
209  }
210  double ss = dot(x, y);
211  return vcl_acos(ss / vcl_sqrt(den));
212  }
213 
214  // Deep copy operator
216  {
217  return new Self(*this);
218  }
219 
220 protected:
221  SAMKernelFunctor(const Self& copy)
222  : Superclass(copy)
223  {
224  *this = copy;
225  }
226 
228  {
230  return *this;
231  }
232 };
233 
238 {
239 public:
242 
244  {
245  this->SetName("RadialSAM");
246  }
248 
249  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
250  {
251  return vcl_exp(-param.gamma * m_Sam(x, y, param));
252  }
253 
254  // Deep copy operator
256  {
257  return new Self(*this);
258  }
259 
260 protected:
262  : Superclass(copy)
263  {
264  *this = copy;
265  }
266 
268  {
270  return *this;
271  }
272 
273 private:
275 };
276 
281 {
282 public:
285 
287  {
288  this->SetName("InverseCosSAM");
289  }
291 
292  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
293  {
294  return 1.0 - vcl_cos(m_Sam(x, y, param));
295  }
296 
297  // Deep copy operator
299  {
300  return new Self(*this);
301  }
302 
303 protected:
305  : Superclass(copy)
306  {
307  *this = copy;
308  }
309 
311  {
313  return *this;
314  }
315 
316 private:
318 };
319 
324 {
325 public:
328 
330  {
331  this->SetName("invMultiQuadraticSAM");
332  this->SetValue<double>("const_coef", 1);
333  }
335 
336  // Deep copy operator
338  {
339  return new Self(*this);
340  }
341 
342  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
343  {
344  double mq = this->GetValue<double>("const_coef") + m_Sam(x, y, param);
345 
346  if (mq == 0.)
347  {
349  }
350  return 1. / sqrt(mq);
351  }
352 
353 protected:
355  : Superclass(copy)
356  {
357  *this = copy;
358  }
359 
361  {
363  return *this;
364  }
365 
366 private:
368 };
369 
374 {
375 public:
378 
380  {
381  this->SetName("KModSAM");
382  this->SetValue<double>("const_coef", 1);
383  }
385 
386  // Deep copy operator
388  {
389  return new Self(*this);
390  }
391 
392  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
393  {
394  double mq = this->GetValue<double>("const_coef") + m_Sam(x, y, param);
395 
396  if (mq == 0.)
397  {
399  }
400  return vcl_exp(param.gamma / mq) - 1.0;
401  }
402 
403 protected:
405  : Superclass(copy)
406  {
407  *this = copy;
408  }
409 
411  {
413  return *this;
414  }
415 
416 private:
418 };
419 
424 {
425 public:
428 
430  {
431  this->SetName("RBF");
432  this->SetValue<double>("gamma_coef", 0.5);
433  }
434  virtual ~RBFKernelFunctor() {}
435 
436  // Deep copy operator
438  {
439  return new Self(*this);
440  }
441 
442  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
443  {
444  double res = this->GetValue<double>("gamma_coef") * m_Custom(x, y, param);
445  return vcl_exp(-res);
446  }
447 
448  virtual double derivative(const svm_node *x,
449  const svm_node *y,
450  const svm_parameter& param,
451  int degree,
452  int index,
453  bool isAtEnd,
454  double constValue) const
455  {
456  double gamma = this->GetValue<double>("gamma_coef");
457  double kernelValue = 0.;
458  double xval = 0.;
459  double yval = 0.;
460  int compt = 0;
461 
462  const svm_node *xtemp = x;
463  const svm_node *ytemp = y;
464 
465  bool stop = false;
466  while (xtemp->index != -1 && ytemp->index != -1 && stop == false)
467  {
468  if (xtemp->index == ytemp->index)
469  {
470  if (compt == index)
471  {
472  xval = xtemp->value;
473  yval = ytemp->value;
474  stop = true;
475  }
476  else
477  {
478  compt++;
479  ++ytemp;
480  ++xtemp;
481  }
482  }
483  else
484  {
485  if (xtemp->index > ytemp->index) ++ytemp;
486  else ++xtemp;
487  }
488  }
489 
490  if (isAtEnd == true)
491  {
492  kernelValue = this->operator ()(x, y, param);
493  }
494  else
495  {
496  kernelValue = constValue;
497  }
498 
499  if (degree < 0)
500  {
501  return 0;
502  }
503  switch (degree)
504  {
505  case 0:
506  return kernelValue;
507  break;
508  case 1:
509  return (-2 * gamma * (yval - xval) * kernelValue);
510  break;
511  default:
512  return (-2 * gamma *
513  ((degree - 1) *
514  this->derivative(x, y, param, degree - 2, index, isAtEnd, constValue) +
515  (yval - xval) * derivative(x, y, param, degree - 1, index, isAtEnd, constValue)));
516  break;
517  }
518 
519  }
520 
521 protected:
522  RBFKernelFunctor(const Self& copy)
523  : Superclass(copy)
524  {
525  *this = copy;
526  }
527 
529  {
531  return *this;
532  }
533 
534 private:
536 };
537 
542 {
543 public:
546 
548  {
549  this->SetName("RBFRBFSAM");
550  this->SetValue<double>("lin_coef", 1.);
551  }
553 
554  // Deep copy operator
556  {
557  return new Self(*this);
558  }
559 
560  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
561  {
562  return (this->GetValue<double>("lin_coef") * vcl_exp(-param.gamma * m_Custom(x, y, param))
563  + (1.0 - this->GetValue<double>("lin_coef")) * m_RadialSam(x, y, param));
564  }
565 
566 protected:
568  : Superclass(copy)
569  {
570  *this = copy;
571  }
572 
574  {
576  return *this;
577  }
578 
579 private:
582 };
583 
588 {
589 public:
592 
594  {
595  this->SetName("polyRBFSAM");
596  this->SetValue<double>("const_coef", 1);
597  this->SetValue<double>("lin_coef", 1);
598  }
600 
601  // Deep copy operator
603  {
604  return new Self(*this);
605  }
606 
607  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
608  {
609  return this->GetValue<double>("const_lin") * vcl_pow(dot(x, y) + this->GetValue<double>("const_coef"), param.degree)
610  + (1.0 - this->GetValue<double>("const_coef"))
611  * m_RadialSam(x, y, param);
612  }
613 
614 protected:
616  : Superclass(copy)
617  {
618  *this = copy;
619  }
620 
622  {
624  return *this;
625  }
626 
627 private:
629 };
630 
635 {
636 public:
639 
641  {
642  this->SetName("RBFDiff");
643  }
645 
646  // Deep copy operator
648  {
649  return new Self(*this);
650  }
651 
652  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
653  {
654 
655  double total = 0.;
656 
657  while (x->index != -1 && y->index != -1)
658  {
659  if (x->index == y->index)
660  {
661  total = total + vcl_exp(-param.gamma * vcl_abs(x->value - y->value));
662  ++x;
663  ++y;
664  }
665  else
666  {
667  if (x->index < y->index) ++x;
668  else ++y;
669  }
670  }
671  return total;
672  }
673 
674 protected:
676  : Superclass(copy)
677  {
678  *this = copy;
679  }
680 
682  {
684  return *this;
685  }
686 
687 };
688 
693 {
694 public:
697 
699  {
700  this->SetName("customLinear");
701  }
703 
704  // Deep copy operator
706  {
707  return new Self(*this);
708  }
709 
710  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
711  {
712  double value=-param.gamma*m_Custom(x, y, param);
713  if(value<-30.0)
714  {
715  return 0.0;
716  }
717 
718  return (vcl_exp(value));
719  }
720 
721 protected:
723  : Superclass(copy)
724  {
725  *this = copy;
726  }
727 
729  {
731  return *this;
732  }
733 
734 private:
736 };
737 
742 {
743 public:
746 
748  {
749  this->SetName("groupedRBF");
750  }
752 
753  // Deep copy operator
755  {
756  return new Self(*this);
757  }
758 
759  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
760  {
761  const char * parameters = param.custom;
762  const char * twoPoints = ":";
763  const char * position = strpbrk (parameters, twoPoints);
764  int twoPointsPosition = 0;
765  double total = 0.;
766  double value = 0.;
767  int numberOfGroups = atoi(parameters);
768 
769  twoPointsPosition = position - parameters;
770 
771  parameters = parameters + twoPointsPosition + 1;
772 
773  int i, j;
774  std::vector<int> begin;
775  begin.resize(numberOfGroups);
776  std::vector<int> end;
777  end.resize(numberOfGroups);
778 
779  for (i = 0; i < numberOfGroups; ++i)
780  {
781  begin[i] = atoi(parameters);
782  position = strpbrk (parameters, twoPoints);
783  twoPointsPosition = position - parameters;
784  parameters = parameters + twoPointsPosition + 1;
785 
786  end[i] = atoi(parameters);
787  parameters = parameters + twoPointsPosition + 1;
788  }
789 
790  const svm_node * xBuff = x;
791  const svm_node * yBuff = y;
792 
793  int sizeX = 0;
794  int sizeY = 0;
795 
796  while (xBuff->index != 1)
797  {
798  ++sizeX;
799  ++xBuff;
800  }
801 
802  while (yBuff->index != 1)
803  {
804  ++sizeY;
805  ++yBuff;
806  }
807 
808  const svm_node * xTemp = NULL;
809  const svm_node * yTemp = NULL;
810 
811  int index = 0;
812 
813  if (sizeX && sizeY)
814  {
815  svm_node* xGroup = new svm_node[sizeX];
816  svm_node* yGroup = new svm_node[sizeY];
817  for (j = 0; j < numberOfGroups; ++j)
818  {
819  xTemp = x;
820  yTemp = y;
821  index = 0;
822 
823  while (xTemp->index != 1 && yTemp->index != 1)
824  {
825  xGroup[index].index = xTemp->index;
826  yGroup[index].index = yTemp->index;
827  if ((xTemp->index < begin[j]) || (xTemp->index > end[j]))
828  {
829  xGroup[index].value = 0;
830  yGroup[index].value = 0;
831  }
832  else
833  {
834  xGroup[index].value = xTemp->value;
835  yGroup[index].value = yTemp->value;
836  }
837  ++index;
838  ++xTemp;
839  ++yTemp;
840  }
841 
842  // value can have different value according to j
843  value = m_CustomLinear(xGroup, yGroup, param);
844 
845  total += value;
846  }
847  delete[] xGroup;
848  delete[] yGroup;
849  }
850 
851  else if ((sizeX > 0) && (sizeY == 0))
852  {
853  svm_node* xGroup = new svm_node[sizeX];
854  svm_node* yGroup = new svm_node[sizeY];
855 
856  for (j = 0; j < numberOfGroups; ++j)
857  {
858  xTemp = x;
859  index = 0;
860 
861  while (xTemp->index != -1)
862  {
863  xGroup[index].index = xTemp->index;
864  yGroup[index].index = xTemp->index;
865 
866  if ((xTemp->index < begin[j]) || (xTemp->index > end[j]))
867  {
868  xGroup[index].value = 0;
869  yGroup[index].value = 0;
870  }
871  else
872  {
873  xGroup[index].value = xTemp->value;
874  yGroup[index].value = 0;
875  }
876  ++index;
877  ++xTemp;
878  }
879 
880  // value can have different value according to j
881  value = m_CustomLinear(xGroup, yGroup, param);
882 
883  total += value;
884  }
885  delete[] xGroup;
886  delete[] yGroup;
887  }
888 
889  else
890  {
891  total = static_cast<double>(numberOfGroups) * m_CustomLinear(x, y, param);
892  }
893 
894  if (xTemp != NULL) delete xTemp;
895  if (yTemp != NULL) delete yTemp;
896 
897  return total;
898  }
899 
900 protected:
902  : Superclass(copy)
903  {
904  *this = copy;
905  }
906 
908  {
910  return *this;
911  }
912 
913 private:
915 };
916 
921 {
922 public:
925 
927  {
928  this->SetName("groupingAdaptive");
929  this->SetValue<double>("lin_coef", 1.);
930  this->SetValue<double>("const_coef", 1.);
931  }
933 
934  // Deep copy operator
936  {
937  return new Self(*this);
938  }
939 
940  virtual double operator ()(const svm_node *x, const svm_node *y, const svm_parameter& param) const
941  {
942  const char * parameters = param.custom;
943  const char * twoPoints = ":";
944  const char * position = strpbrk (parameters, twoPoints);
945  int twoPointsPosition = 0;
946  double total = 0.;
947  double value = 0.;
948  int numberOfGroups = atoi(parameters);
949 
950  twoPointsPosition = position - parameters;
951 
952  parameters = parameters + twoPointsPosition + 1;
953 
954  int i, j;
955 
956  std::vector<int> begin;
957  begin.resize(numberOfGroups);
958  std::vector<int> end;
959  end.resize(numberOfGroups);
960 
961  for (i = 0; i < numberOfGroups; ++i)
962  {
963  begin[i] = atoi(parameters);
964  position = strpbrk (parameters, twoPoints);
965  twoPointsPosition = position - parameters;
966  parameters = parameters + twoPointsPosition + 1;
967 
968  end[i] = atoi(parameters);
969  parameters = parameters + twoPointsPosition + 1;
970  }
971 
972  const svm_node * xBuff = x;
973  const svm_node * yBuff = y;
974 
975  int sizeX = 0;
976  int sizeY = 0;
977 
978  while (xBuff->index != 1)
979  {
980  ++sizeX;
981  ++xBuff;
982  }
983 
984  while (yBuff->index != 1)
985  {
986  ++sizeY;
987  ++yBuff;
988  }
989 
990  const svm_node * xTemp = NULL;
991  const svm_node * yTemp = NULL;
992 
993  int index = 0;
994 
995  if (sizeX && sizeY)
996  {
997  svm_node* xGroup = new svm_node[sizeX];
998  svm_node* yGroup = new svm_node[sizeY];
999 
1000  for (j = 0; j < numberOfGroups; ++j)
1001  {
1002  xTemp = x;
1003  yTemp = y;
1004  index = 0;
1005 
1006  while (xTemp->index != 1 && yTemp->index != 1)
1007  {
1008  xGroup[index].index = xTemp->index;
1009  yGroup[index].index = yTemp->index;
1010  if ((xTemp->index < begin[j]) || (xTemp->index > end[j]))
1011  {
1012  xGroup[index].value = 0;
1013  yGroup[index].value = 0;
1014  }
1015  else
1016  {
1017  xGroup[index].value = xTemp->value;
1018  yGroup[index].value = yTemp->value;
1019  }
1020  ++index;
1021  ++xTemp;
1022  ++yTemp;
1023  }
1024 
1025  // value can have different value according to j
1026  value = vcl_pow(this->GetValue<double>("lin_coef") * dot(xGroup, yGroup) + this->GetValue<double>(
1027  "const_coef"), static_cast<double>(param.degree));
1028 
1029  total += value;
1030  }
1031  delete[] xGroup;
1032  delete[] yGroup;
1033  }
1034 
1035  else if ((sizeX > 0) && (sizeY == 0))
1036  {
1037  svm_node* xGroup = new svm_node[sizeX];
1038  svm_node* yGroup = new svm_node[sizeY];
1039 
1040  for (j = 0; j < numberOfGroups; ++j)
1041  {
1042  xTemp = x;
1043  index = 0;
1044 
1045  while (xTemp->index != 1)
1046  {
1047  xGroup[index].index = xTemp->index;
1048  yGroup[index].index = xTemp->index;
1049 
1050  if ((xTemp->index < begin[j]) || (xTemp->index > end[j]))
1051  {
1052  xGroup[index].value = 0;
1053  yGroup[index].value = 0;
1054  }
1055  else
1056  {
1057  xGroup[index].value = xTemp->value;
1058  yGroup[index].value = 0;
1059  }
1060  ++index;
1061  ++xTemp;
1062  }
1063 
1064  // value can have different value according to j
1065  value = vcl_pow(this->GetValue<double>("lin_coef") * dot(xGroup, yGroup) + this->GetValue<double>(
1066  "const_coef"), static_cast<double>(param.degree));
1067 
1068  total += value;
1069  }
1070  delete[] xGroup;
1071  delete[] yGroup;
1072  }
1073 
1074  else
1075  {
1076  CustomKernelFunctor custom;
1077  total = static_cast<double>(numberOfGroups) * custom(x, y, param);
1078  }
1079 
1080  return total;
1081  }
1082 
1083 
1084 protected:
1086  : Superclass(copy)
1087  {
1088  *this = copy;
1089  }
1090 
1092  {
1093  Superclass::operator =(copy);
1094  return *this;
1095  }
1096 
1097 };
1098 
1099 }
1100 #endif

Generated at Sat Mar 8 2014 16:07:40 for Orfeo Toolbox with doxygen 1.8.3.1