kullbackleiblerdistance.h

/* 
 * Copyright (C) 2006-2013  Music Technology Group - Universitat Pompeu Fabra
 *
 * This file is part of Gaia
 * 
 * Gaia is free software: you can redistribute it and/or modify it under 
 * the terms of the GNU Affero General Public License as published by the Free 
 * Software Foundation (FSF), either version 3 of the License, or (at your 
 * option) any later version.
 * 
 * This program is distributed in the hope that it will be useful, but WITHOUT 
 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 
 * FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more 
 * details.
 * 
 * You should have received a copy of the Affero GNU General Public License
 * version 3 along with this program.  If not, see http://www.gnu.org/licenses/
 */

#ifndef GAIA_KULLBACKLEIBLERDISTANCE_H
#define GAIA_KULLBACKLEIBLERDISTANCE_H

#include "distancefunction.h"

namespace gaia2 {

class KullbackLeiblerDistance : public DistanceFunction {
 public:
  KullbackLeiblerDistance(const PointLayout& layout, const ParameterMap& params);
  Real operator()(const Point& p1, const Point& p2, int seg1, int seg2) const;

 protected:
  int _meanIdx, _covIdx, _icovIdx;
  int _meanIdxEnd, _covIdxEnd, _icovIdxEnd;
};


inline Real traceMatrixProduct(const Real* m1, const Real* m2, int size) {
  Real result = 0.0;
  for (int i=0; i<size; i++) {
    for (int j=0; j<size; j++) {
      result += m1[i*size + j] * m2[j*size + i];
    }
  }
  return result;
}


inline Real traceMatrixProductMean(const Real* m1, const Real* m2,
                                   const Real* v1, const Real* v2,
                                   int size) {
  Real result = 0.0;

  // we want this to be on the stack to avoid mem allocation, this function
  // needs to be real fast
  Q_ASSERT(size < 64);
  Real v[64];

  for (int i=0; i<size; i++) v[i] = v1[i] - v2[i];

  for (int i=0; i<size; i++) {
    for (int j=0; j<size; j++) {
      result += (m1[i*size + j] + m2[i*size + j]) * v[j] * v[i];
    }
  }

  return result;
}



} // namespace gaia2

#endif // GAIA_KULLBACKLEIBLERDISTANCE_H