blob: 452c17f4e82c8d5f184564884030a3a2cbd6c965 (
plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
|
// Copyright (c) the JPEG XL Project Authors. All rights reserved.
//
// Use of this source code is governed by a BSD-style
// license that can be found in the LICENSE file.
#include "lib/jxl/enc_linalg.h"
#include <cmath>
#include "lib/jxl/base/status.h"
namespace jxl {
void ConvertToDiagonal(const Matrix2x2& A, Vector2& diag, Matrix2x2& U) {
#if JXL_ENABLE_ASSERT
// Check A is symmetric.
JXL_ASSERT(std::abs(A[0][1] - A[1][0]) < 1e-15);
#endif
if (std::abs(A[0][1]) < 1e-15) {
// Already diagonal.
diag[0] = A[0][0];
diag[1] = A[1][1];
U[0][0] = U[1][1] = 1.0;
U[0][1] = U[1][0] = 0.0;
return;
}
double b = -(A[0][0] + A[1][1]);
double c = A[0][0] * A[1][1] - A[0][1] * A[0][1];
double d = b * b - 4.0 * c;
double sqd = std::sqrt(d);
double l1 = (-b - sqd) * 0.5;
double l2 = (-b + sqd) * 0.5;
Vector2 v1 = {A[0][0] - l1, A[1][0]};
double v1n = 1.0 / std::hypot(v1[0], v1[1]);
v1[0] = v1[0] * v1n;
v1[1] = v1[1] * v1n;
diag[0] = l1;
diag[1] = l2;
U[0][0] = v1[1];
U[0][1] = -v1[0];
U[1][0] = v1[0];
U[1][1] = v1[1];
}
} // namespace jxl
|
