00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00028 #include "avassert.h"
00029
00030 #include <limits.h>
00031
00032 #include "common.h"
00033 #include "mathematics.h"
00034 #include "rational.h"
00035
00036 int av_reduce(int *dst_num, int *dst_den,
00037 int64_t num, int64_t den, int64_t max)
00038 {
00039 AVRational a0 = { 0, 1 }, a1 = { 1, 0 };
00040 int sign = (num < 0) ^ (den < 0);
00041 int64_t gcd = av_gcd(FFABS(num), FFABS(den));
00042
00043 if (gcd) {
00044 num = FFABS(num) / gcd;
00045 den = FFABS(den) / gcd;
00046 }
00047 if (num <= max && den <= max) {
00048 a1 = (AVRational) { num, den };
00049 den = 0;
00050 }
00051
00052 while (den) {
00053 uint64_t x = num / den;
00054 int64_t next_den = num - den * x;
00055 int64_t a2n = x * a1.num + a0.num;
00056 int64_t a2d = x * a1.den + a0.den;
00057
00058 if (a2n > max || a2d > max) {
00059 if (a1.num) x = (max - a0.num) / a1.num;
00060 if (a1.den) x = FFMIN(x, (max - a0.den) / a1.den);
00061
00062 if (den * (2 * x * a1.den + a0.den) > num * a1.den)
00063 a1 = (AVRational) { x * a1.num + a0.num, x * a1.den + a0.den };
00064 break;
00065 }
00066
00067 a0 = a1;
00068 a1 = (AVRational) { a2n, a2d };
00069 num = den;
00070 den = next_den;
00071 }
00072 av_assert2(av_gcd(a1.num, a1.den) <= 1U);
00073
00074 *dst_num = sign ? -a1.num : a1.num;
00075 *dst_den = a1.den;
00076
00077 return den == 0;
00078 }
00079
00080 AVRational av_mul_q(AVRational b, AVRational c)
00081 {
00082 av_reduce(&b.num, &b.den,
00083 b.num * (int64_t) c.num,
00084 b.den * (int64_t) c.den, INT_MAX);
00085 return b;
00086 }
00087
00088 AVRational av_div_q(AVRational b, AVRational c)
00089 {
00090 return av_mul_q(b, (AVRational) { c.den, c.num });
00091 }
00092
00093 AVRational av_add_q(AVRational b, AVRational c) {
00094 av_reduce(&b.num, &b.den,
00095 b.num * (int64_t) c.den +
00096 c.num * (int64_t) b.den,
00097 b.den * (int64_t) c.den, INT_MAX);
00098 return b;
00099 }
00100
00101 AVRational av_sub_q(AVRational b, AVRational c)
00102 {
00103 return av_add_q(b, (AVRational) { -c.num, c.den });
00104 }
00105
00106 AVRational av_d2q(double d, int max)
00107 {
00108 AVRational a;
00109 #define LOG2 0.69314718055994530941723212145817656807550013436025
00110 int exponent;
00111 int64_t den;
00112 if (isnan(d))
00113 return (AVRational) { 0,0 };
00114 if (isinf(d))
00115 return (AVRational) { d < 0 ? -1 : 1, 0 };
00116 exponent = FFMAX( (int)(log(fabs(d) + 1e-20)/LOG2), 0);
00117 den = 1LL << (61 - exponent);
00118 av_reduce(&a.num, &a.den, (int64_t)(d * den + 0.5), den, max);
00119
00120 return a;
00121 }
00122
00123 int av_nearer_q(AVRational q, AVRational q1, AVRational q2)
00124 {
00125
00126 int64_t a = q1.num * (int64_t)q2.den + q2.num * (int64_t)q1.den;
00127 int64_t b = 2 * (int64_t)q1.den * q2.den;
00128
00129
00130 int64_t x_up = av_rescale_rnd(a, q.den, b, AV_ROUND_UP);
00131
00132
00133 int64_t x_down = av_rescale_rnd(a, q.den, b, AV_ROUND_DOWN);
00134
00135 return ((x_up > q.num) - (x_down < q.num)) * av_cmp_q(q2, q1);
00136 }
00137
00138 int av_find_nearest_q_idx(AVRational q, const AVRational* q_list)
00139 {
00140 int i, nearest_q_idx = 0;
00141 for (i = 0; q_list[i].den; i++)
00142 if (av_nearer_q(q, q_list[i], q_list[nearest_q_idx]) > 0)
00143 nearest_q_idx = i;
00144
00145 return nearest_q_idx;
00146 }
00147
00148 #ifdef TEST
00149 int main(void)
00150 {
00151 AVRational a,b;
00152 for (a.num = -2; a.num <= 2; a.num++) {
00153 for (a.den = -2; a.den <= 2; a.den++) {
00154 for (b.num = -2; b.num <= 2; b.num++) {
00155 for (b.den = -2; b.den <= 2; b.den++) {
00156 int c = av_cmp_q(a,b);
00157 double d = av_q2d(a) == av_q2d(b) ?
00158 0 : (av_q2d(a) - av_q2d(b));
00159 if (d > 0) d = 1;
00160 else if (d < 0) d = -1;
00161 else if (d != d) d = INT_MIN;
00162 if (c != d)
00163 av_log(0, AV_LOG_ERROR, "%d/%d %d/%d, %d %f\n", a.num,
00164 a.den, b.num, b.den, c,d);
00165 }
00166 }
00167 }
00168 }
00169 return 0;
00170 }
00171 #endif