2006-07-14 13:03:09 +03:00
|
|
|
/*
|
|
|
|
* linear least squares model
|
|
|
|
*
|
|
|
|
* Copyright (c) 2006 Michael Niedermayer <michaelni@gmx.at>
|
|
|
|
*
|
2006-10-07 18:30:46 +03:00
|
|
|
* This file is part of FFmpeg.
|
|
|
|
*
|
|
|
|
* FFmpeg is free software; you can redistribute it and/or
|
2006-07-14 13:03:09 +03:00
|
|
|
* modify it under the terms of the GNU Lesser General Public
|
|
|
|
* License as published by the Free Software Foundation; either
|
2006-10-07 18:30:46 +03:00
|
|
|
* version 2.1 of the License, or (at your option) any later version.
|
2006-07-14 13:03:09 +03:00
|
|
|
*
|
2006-10-07 18:30:46 +03:00
|
|
|
* FFmpeg is distributed in the hope that it will be useful,
|
2006-07-14 13:03:09 +03:00
|
|
|
* but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
|
|
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
|
|
|
|
* Lesser General Public License for more details.
|
|
|
|
*
|
|
|
|
* You should have received a copy of the GNU Lesser General Public
|
2006-10-07 18:30:46 +03:00
|
|
|
* License along with FFmpeg; if not, write to the Free Software
|
2007-07-05 13:40:25 +03:00
|
|
|
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
|
2006-07-14 13:03:09 +03:00
|
|
|
*/
|
|
|
|
|
|
|
|
/**
|
2010-04-20 17:45:34 +03:00
|
|
|
* @file
|
2006-07-14 13:03:09 +03:00
|
|
|
* linear least squares model
|
|
|
|
*/
|
|
|
|
|
|
|
|
#include <math.h>
|
|
|
|
#include <string.h>
|
|
|
|
|
2013-04-30 19:31:25 +03:00
|
|
|
#include "attributes.h"
|
2015-12-10 22:49:30 +02:00
|
|
|
#include "internal.h"
|
2013-02-28 23:31:26 +03:00
|
|
|
#include "version.h"
|
2014-08-10 00:03:26 +03:00
|
|
|
#include "lls.h"
|
2006-07-14 13:03:09 +03:00
|
|
|
|
2014-09-28 16:44:38 +03:00
|
|
|
static void update_lls(LLSModel *m, const double *var)
|
2011-07-08 14:52:12 +03:00
|
|
|
{
|
|
|
|
int i, j;
|
2006-07-14 13:03:09 +03:00
|
|
|
|
2011-07-08 14:52:12 +03:00
|
|
|
for (i = 0; i <= m->indep_count; i++) {
|
|
|
|
for (j = i; j <= m->indep_count; j++) {
|
|
|
|
m->covariance[i][j] += var[i] * var[j];
|
2006-07-14 13:03:09 +03:00
|
|
|
}
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
2014-08-10 00:03:26 +03:00
|
|
|
void avpriv_solve_lls(LLSModel *m, double threshold, unsigned short min_order)
|
2011-07-08 14:52:12 +03:00
|
|
|
{
|
|
|
|
int i, j, k;
|
2013-06-19 00:30:42 +03:00
|
|
|
double (*factor)[MAX_VARS_ALIGN] = (void *) &m->covariance[1][0];
|
|
|
|
double (*covar) [MAX_VARS_ALIGN] = (void *) &m->covariance[1][1];
|
2011-07-08 14:52:12 +03:00
|
|
|
double *covar_y = m->covariance[0];
|
|
|
|
int count = m->indep_count;
|
|
|
|
|
|
|
|
for (i = 0; i < count; i++) {
|
|
|
|
for (j = i; j < count; j++) {
|
|
|
|
double sum = covar[i][j];
|
|
|
|
|
avutil/lls: speed up performance of solve_lls
This is a trivial rewrite of the loops that results in better
prefetching and associated cache efficiency. Essentially, the problem is
that modern prefetching logic is based on finite state Markov memory, a reasonable
assumption that is used elsewhere in CPU's in for instance branch
predictors.
Surrounding loops all iterate forward through the array, making the
predictor think of prefetching in the forward direction, but the
intermediate loop is unnecessarily in the backward direction.
Speedup is nontrivial. Benchmarks obtained by 10^6 iterations within
solve_lls, with START/STOP_TIMER. File is tests/data/fate/flac-16-lpc-cholesky.err.
Hardware: x86-64, Haswell, GNU/Linux.
new:
17291 decicycles in solve_lls, 2096706 runs, 446 skips
17255 decicycles in solve_lls, 4193657 runs, 647 skips
17231 decicycles in solve_lls, 8384997 runs, 3611 skips
17189 decicycles in solve_lls,16771010 runs, 6206 skips
17132 decicycles in solve_lls,33544757 runs, 9675 skips
17092 decicycles in solve_lls,67092404 runs, 16460 skips
17058 decicycles in solve_lls,134188213 runs, 29515 skips
old:
18009 decicycles in solve_lls, 2096665 runs, 487 skips
17805 decicycles in solve_lls, 4193320 runs, 984 skips
17779 decicycles in solve_lls, 8386855 runs, 1753 skips
18289 decicycles in solve_lls,16774280 runs, 2936 skips
18158 decicycles in solve_lls,33548104 runs, 6328 skips
18420 decicycles in solve_lls,67091793 runs, 17071 skips
18310 decicycles in solve_lls,134187219 runs, 30509 skips
Reviewed-by: Michael Niedermayer <michael@niedermayer.cc>
Signed-off-by: Ganesh Ajjanagadde <gajjanagadde@gmail.com>
2015-11-24 02:05:54 +02:00
|
|
|
for (k = 0; k <= i-1; k++)
|
2011-07-08 14:52:12 +03:00
|
|
|
sum -= factor[i][k] * factor[j][k];
|
|
|
|
|
|
|
|
if (i == j) {
|
|
|
|
if (sum < threshold)
|
|
|
|
sum = 1.0;
|
|
|
|
factor[i][i] = sqrt(sum);
|
|
|
|
} else {
|
|
|
|
factor[j][i] = sum / factor[i][i];
|
|
|
|
}
|
2006-07-14 13:03:09 +03:00
|
|
|
}
|
|
|
|
}
|
2011-07-08 14:52:12 +03:00
|
|
|
|
|
|
|
for (i = 0; i < count; i++) {
|
|
|
|
double sum = covar_y[i + 1];
|
|
|
|
|
avutil/lls: speed up performance of solve_lls
This is a trivial rewrite of the loops that results in better
prefetching and associated cache efficiency. Essentially, the problem is
that modern prefetching logic is based on finite state Markov memory, a reasonable
assumption that is used elsewhere in CPU's in for instance branch
predictors.
Surrounding loops all iterate forward through the array, making the
predictor think of prefetching in the forward direction, but the
intermediate loop is unnecessarily in the backward direction.
Speedup is nontrivial. Benchmarks obtained by 10^6 iterations within
solve_lls, with START/STOP_TIMER. File is tests/data/fate/flac-16-lpc-cholesky.err.
Hardware: x86-64, Haswell, GNU/Linux.
new:
17291 decicycles in solve_lls, 2096706 runs, 446 skips
17255 decicycles in solve_lls, 4193657 runs, 647 skips
17231 decicycles in solve_lls, 8384997 runs, 3611 skips
17189 decicycles in solve_lls,16771010 runs, 6206 skips
17132 decicycles in solve_lls,33544757 runs, 9675 skips
17092 decicycles in solve_lls,67092404 runs, 16460 skips
17058 decicycles in solve_lls,134188213 runs, 29515 skips
old:
18009 decicycles in solve_lls, 2096665 runs, 487 skips
17805 decicycles in solve_lls, 4193320 runs, 984 skips
17779 decicycles in solve_lls, 8386855 runs, 1753 skips
18289 decicycles in solve_lls,16774280 runs, 2936 skips
18158 decicycles in solve_lls,33548104 runs, 6328 skips
18420 decicycles in solve_lls,67091793 runs, 17071 skips
18310 decicycles in solve_lls,134187219 runs, 30509 skips
Reviewed-by: Michael Niedermayer <michael@niedermayer.cc>
Signed-off-by: Ganesh Ajjanagadde <gajjanagadde@gmail.com>
2015-11-24 02:05:54 +02:00
|
|
|
for (k = 0; k <= i-1; k++)
|
2011-07-08 14:52:12 +03:00
|
|
|
sum -= factor[i][k] * m->coeff[0][k];
|
|
|
|
|
|
|
|
m->coeff[0][i] = sum / factor[i][i];
|
2006-07-14 13:03:09 +03:00
|
|
|
}
|
|
|
|
|
2011-07-08 14:52:12 +03:00
|
|
|
for (j = count - 1; j >= min_order; j--) {
|
|
|
|
for (i = j; i >= 0; i--) {
|
|
|
|
double sum = m->coeff[0][i];
|
|
|
|
|
|
|
|
for (k = i + 1; k <= j; k++)
|
|
|
|
sum -= factor[k][i] * m->coeff[j][k];
|
|
|
|
|
|
|
|
m->coeff[j][i] = sum / factor[i][i];
|
2006-07-16 02:43:38 +03:00
|
|
|
}
|
2006-07-14 13:03:09 +03:00
|
|
|
|
2011-07-08 14:52:12 +03:00
|
|
|
m->variance[j] = covar_y[0];
|
|
|
|
|
|
|
|
for (i = 0; i <= j; i++) {
|
|
|
|
double sum = m->coeff[j][i] * covar[i][i] - 2 * covar_y[i + 1];
|
|
|
|
|
|
|
|
for (k = 0; k < i; k++)
|
|
|
|
sum += 2 * m->coeff[j][k] * covar[k][i];
|
|
|
|
|
|
|
|
m->variance[j] += m->coeff[j][i] * sum;
|
2006-07-16 02:43:38 +03:00
|
|
|
}
|
2006-07-14 13:03:09 +03:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
2014-09-28 16:44:38 +03:00
|
|
|
static double evaluate_lls(LLSModel *m, const double *param, int order)
|
2011-07-08 14:52:12 +03:00
|
|
|
{
|
2006-07-14 13:03:09 +03:00
|
|
|
int i;
|
2011-07-08 14:52:12 +03:00
|
|
|
double out = 0;
|
2006-07-14 13:03:09 +03:00
|
|
|
|
2011-07-08 14:52:12 +03:00
|
|
|
for (i = 0; i <= order; i++)
|
|
|
|
out += param[i] * m->coeff[order][i];
|
2006-07-14 13:03:09 +03:00
|
|
|
|
|
|
|
return out;
|
|
|
|
}
|
|
|
|
|
2014-08-10 00:03:26 +03:00
|
|
|
av_cold void avpriv_init_lls(LLSModel *m, int indep_count)
|
2013-06-19 00:30:41 +03:00
|
|
|
{
|
2014-08-10 00:03:26 +03:00
|
|
|
memset(m, 0, sizeof(LLSModel));
|
2013-06-19 00:30:41 +03:00
|
|
|
m->indep_count = indep_count;
|
|
|
|
m->update_lls = update_lls;
|
|
|
|
m->evaluate_lls = evaluate_lls;
|
2013-06-19 00:30:42 +03:00
|
|
|
if (ARCH_X86)
|
|
|
|
ff_init_lls_x86(m);
|
2013-06-19 00:30:41 +03:00
|
|
|
}
|
|
|
|
|
2006-07-14 13:03:09 +03:00
|
|
|
#ifdef TEST
|
|
|
|
|
|
|
|
#include <stdio.h>
|
2011-07-08 17:01:05 +03:00
|
|
|
#include <limits.h>
|
|
|
|
#include "lfg.h"
|
2006-07-14 13:03:09 +03:00
|
|
|
|
2011-07-08 14:52:12 +03:00
|
|
|
int main(void)
|
|
|
|
{
|
2014-08-10 00:03:26 +03:00
|
|
|
LLSModel m;
|
2006-07-16 02:43:38 +03:00
|
|
|
int i, order;
|
2011-07-08 17:01:05 +03:00
|
|
|
AVLFG lfg;
|
2006-07-14 13:03:09 +03:00
|
|
|
|
2011-07-08 17:01:05 +03:00
|
|
|
av_lfg_init(&lfg, 1);
|
2014-08-10 00:03:26 +03:00
|
|
|
avpriv_init_lls(&m, 3);
|
2006-07-14 13:03:09 +03:00
|
|
|
|
2011-07-08 14:52:12 +03:00
|
|
|
for (i = 0; i < 100; i++) {
|
2013-06-19 00:30:42 +03:00
|
|
|
LOCAL_ALIGNED(32, double, var, [4]);
|
2008-01-09 00:48:26 +02:00
|
|
|
double eval;
|
2011-07-08 14:52:12 +03:00
|
|
|
|
2016-03-21 16:19:30 +02:00
|
|
|
var[0] = (av_lfg_get(&lfg) / (double) UINT_MAX - 0.5) * 2;
|
2011-07-08 17:01:05 +03:00
|
|
|
var[1] = var[0] + av_lfg_get(&lfg) / (double) UINT_MAX - 0.5;
|
|
|
|
var[2] = var[1] + av_lfg_get(&lfg) / (double) UINT_MAX - 0.5;
|
|
|
|
var[3] = var[2] + av_lfg_get(&lfg) / (double) UINT_MAX - 0.5;
|
2013-06-19 00:30:41 +03:00
|
|
|
m.update_lls(&m, var);
|
2014-08-10 00:03:26 +03:00
|
|
|
avpriv_solve_lls(&m, 0.001, 0);
|
2011-07-08 14:52:12 +03:00
|
|
|
for (order = 0; order < 3; order++) {
|
2013-06-19 00:30:41 +03:00
|
|
|
eval = m.evaluate_lls(&m, var + 1, order);
|
2009-03-20 13:40:05 +02:00
|
|
|
printf("real:%9f order:%d pred:%9f var:%f coeffs:%f %9f %9f\n",
|
2011-07-08 14:52:12 +03:00
|
|
|
var[0], order, eval, sqrt(m.variance[order] / (i + 1)),
|
|
|
|
m.coeff[order][0], m.coeff[order][1],
|
|
|
|
m.coeff[order][2]);
|
2006-07-16 02:43:38 +03:00
|
|
|
}
|
2006-07-14 13:03:09 +03:00
|
|
|
}
|
|
|
|
return 0;
|
|
|
|
}
|
|
|
|
|
|
|
|
#endif
|