FFmpeg
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
rational.c
Go to the documentation of this file.
1 /*
2  * rational numbers
3  * Copyright (c) 2003 Michael Niedermayer <michaelni@gmx.at>
4  *
5  * This file is part of FFmpeg.
6  *
7  * FFmpeg is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU Lesser General Public
9  * License as published by the Free Software Foundation; either
10  * version 2.1 of the License, or (at your option) any later version.
11  *
12  * FFmpeg is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15  * Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with FFmpeg; if not, write to the Free Software
19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20  */
21 
22 /**
23  * @file
24  * rational numbers
25  * @author Michael Niedermayer <michaelni@gmx.at>
26  */
27 
28 #include "avassert.h"
29 #include <limits.h>
30 
31 #include "common.h"
32 #include "mathematics.h"
33 #include "rational.h"
34 
35 int av_reduce(int *dst_num, int *dst_den,
36  int64_t num, int64_t den, int64_t max)
37 {
38  AVRational a0 = { 0, 1 }, a1 = { 1, 0 };
39  int sign = (num < 0) ^ (den < 0);
40  int64_t gcd = av_gcd(FFABS(num), FFABS(den));
41 
42  if (gcd) {
43  num = FFABS(num) / gcd;
44  den = FFABS(den) / gcd;
45  }
46  if (num <= max && den <= max) {
47  a1 = (AVRational) { num, den };
48  den = 0;
49  }
50 
51  while (den) {
52  uint64_t x = num / den;
53  int64_t next_den = num - den * x;
54  int64_t a2n = x * a1.num + a0.num;
55  int64_t a2d = x * a1.den + a0.den;
56 
57  if (a2n > max || a2d > max) {
58  if (a1.num) x = (max - a0.num) / a1.num;
59  if (a1.den) x = FFMIN(x, (max - a0.den) / a1.den);
60 
61  if (den * (2 * x * a1.den + a0.den) > num * a1.den)
62  a1 = (AVRational) { x * a1.num + a0.num, x * a1.den + a0.den };
63  break;
64  }
65 
66  a0 = a1;
67  a1 = (AVRational) { a2n, a2d };
68  num = den;
69  den = next_den;
70  }
71  av_assert2(av_gcd(a1.num, a1.den) <= 1U);
72  av_assert2(a1.num <= max && a1.den <= max);
73 
74  *dst_num = sign ? -a1.num : a1.num;
75  *dst_den = a1.den;
76 
77  return den == 0;
78 }
79 
81 {
82  av_reduce(&b.num, &b.den,
83  b.num * (int64_t) c.num,
84  b.den * (int64_t) c.den, INT_MAX);
85  return b;
86 }
87 
89 {
90  return av_mul_q(b, (AVRational) { c.den, c.num });
91 }
92 
94  av_reduce(&b.num, &b.den,
95  b.num * (int64_t) c.den +
96  c.num * (int64_t) b.den,
97  b.den * (int64_t) c.den, INT_MAX);
98  return b;
99 }
100 
102 {
103  return av_add_q(b, (AVRational) { -c.num, c.den });
104 }
105 
106 AVRational av_d2q(double d, int max)
107 {
108  AVRational a;
109  int exponent;
110  int64_t den;
111  if (isnan(d))
112  return (AVRational) { 0,0 };
113  if (fabs(d) > INT_MAX + 3LL)
114  return (AVRational) { d < 0 ? -1 : 1, 0 };
115  frexp(d, &exponent);
116  exponent = FFMAX(exponent-1, 0);
117  den = 1LL << (61 - exponent);
118  // (int64_t)rint() and llrint() do not work with gcc on ia64 and sparc64,
119  // see Ticket2713 for affected gcc/glibc versions
120  av_reduce(&a.num, &a.den, floor(d * den + 0.5), den, max);
121  if ((!a.num || !a.den) && d && max>0 && max<INT_MAX)
122  av_reduce(&a.num, &a.den, floor(d * den + 0.5), den, INT_MAX);
123 
124  return a;
125 }
126 
128 {
129  /* n/d is q, a/b is the median between q1 and q2 */
130  int64_t a = q1.num * (int64_t)q2.den + q2.num * (int64_t)q1.den;
131  int64_t b = 2 * (int64_t)q1.den * q2.den;
132 
133  /* rnd_up(a*d/b) > n => a*d/b > n */
134  int64_t x_up = av_rescale_rnd(a, q.den, b, AV_ROUND_UP);
135 
136  /* rnd_down(a*d/b) < n => a*d/b < n */
137  int64_t x_down = av_rescale_rnd(a, q.den, b, AV_ROUND_DOWN);
138 
139  return ((x_up > q.num) - (x_down < q.num)) * av_cmp_q(q2, q1);
140 }
141 
143 {
144  int i, nearest_q_idx = 0;
145  for (i = 0; q_list[i].den; i++)
146  if (av_nearer_q(q, q_list[i], q_list[nearest_q_idx]) > 0)
147  nearest_q_idx = i;
148 
149  return nearest_q_idx;
150 }
151 
153  int64_t n;
154  int shift;
155  int sign = 0;
156 
157  if (q.den < 0) {
158  q.den *= -1;
159  q.num *= -1;
160  }
161  if (q.num < 0) {
162  q.num *= -1;
163  sign = 1;
164  }
165 
166  if (!q.num && !q.den) return 0xFFC00000;
167  if (!q.num) return 0;
168  if (!q.den) return 0x7F800000 | (q.num & 0x80000000);
169 
170  shift = 23 + av_log2(q.den) - av_log2(q.num);
171  if (shift >= 0) n = av_rescale(q.num, 1LL<<shift, q.den);
172  else n = av_rescale(q.num, 1, ((int64_t)q.den) << -shift);
173 
174  shift -= n >= (1<<24);
175  shift += n < (1<<23);
176 
177  if (shift >= 0) n = av_rescale(q.num, 1LL<<shift, q.den);
178  else n = av_rescale(q.num, 1, ((int64_t)q.den) << -shift);
179 
180  av_assert1(n < (1<<24));
181  av_assert1(n >= (1<<23));
182 
183  return sign<<31 | (150-shift)<<23 | (n - (1<<23));
184 }
static int shift(int a, int b)
Definition: sonic.c:82
AVRational av_div_q(AVRational b, AVRational c)
Divide one rational by another.
Definition: rational.c:88
#define a0
Definition: regdef.h:46
int num
Numerator.
Definition: rational.h:59
const char * b
Definition: vf_curves.c:113
int av_log2(unsigned v)
Definition: intmath.c:26
#define a1
Definition: regdef.h:47
static const uint8_t q1[256]
Definition: twofish.c:96
uint32_t av_q2intfloat(AVRational q)
Convert an AVRational to a IEEE 32-bit float expressed in fixed-point format.
Definition: rational.c:152
AVRational av_sub_q(AVRational b, AVRational c)
Subtract one rational from another.
Definition: rational.c:101
Round toward +infinity.
Definition: mathematics.h:83
#define av_assert2(cond)
assert() equivalent, that does lie in speed critical code.
Definition: avassert.h:64
int av_reduce(int *dst_num, int *dst_den, int64_t num, int64_t den, int64_t max)
Reduce a fraction.
Definition: rational.c:35
#define U(x)
Definition: vp56_arith.h:37
simple assert() macros that are a bit more flexible than ISO C assert().
int64_t av_gcd(int64_t a, int64_t b)
Compute the greatest common divisor of two integer operands.
Definition: mathematics.c:37
#define FFMAX(a, b)
Definition: common.h:94
#define av_assert1(cond)
assert() equivalent, that does not lie in speed critical code.
Definition: avassert.h:53
int64_t av_rescale(int64_t a, int64_t b, int64_t c)
Rescale a 64-bit integer with rounding to nearest.
Definition: mathematics.c:129
#define FFMIN(a, b)
Definition: common.h:96
#define FFABS(a)
Absolute value, Note, INT_MIN / INT64_MIN result in undefined behavior as they are not representable ...
Definition: common.h:72
int n
Definition: avisynth_c.h:684
int64_t av_rescale_rnd(int64_t a, int64_t b, int64_t c, enum AVRounding rnd)
Rescale a 64-bit integer with specified rounding.
Definition: mathematics.c:58
int av_find_nearest_q_idx(AVRational q, const AVRational *q_list)
Find the value in a list of rationals nearest a given reference rational.
Definition: rational.c:142
Rational number (pair of numerator and denominator).
Definition: rational.h:58
#define isnan(x)
Definition: libm.h:340
Round toward -infinity.
Definition: mathematics.h:82
AVRational av_d2q(double d, int max)
Convert a double precision floating point number to a rational.
Definition: rational.c:106
static int av_cmp_q(AVRational a, AVRational b)
Compare two rationals.
Definition: rational.h:89
common internal and external API header
if(ret< 0)
Definition: vf_mcdeint.c:282
Utilties for rational number calculation.
static double c[64]
AVRational av_add_q(AVRational b, AVRational c)
Add two rationals.
Definition: rational.c:93
int av_nearer_q(AVRational q, AVRational q1, AVRational q2)
Find which of the two rationals is closer to another rational.
Definition: rational.c:127
int den
Denominator.
Definition: rational.h:60
AVRational av_mul_q(AVRational b, AVRational c)
Multiply two rationals.
Definition: rational.c:80