1 /*
2 * Copyright (C) 2011 The Android Open Source Project
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16
17 #define __STDC_LIMIT_MACROS
18
19 #include <assert.h>
20 #include <stdint.h>
21
22 #include <utils/LinearTransform.h>
23
24 namespace android {
25
ABS(T x)26 template<class T> static inline T ABS(T x) { return (x < 0) ? -x : x; }
27
28 // Static math methods involving linear transformations
scale_u64_to_u64(uint64_t val,uint32_t N,uint32_t D,uint64_t * res,bool round_up_not_down)29 static bool scale_u64_to_u64(
30 uint64_t val,
31 uint32_t N,
32 uint32_t D,
33 uint64_t* res,
34 bool round_up_not_down) {
35 uint64_t tmp1, tmp2;
36 uint32_t r;
37
38 assert(res);
39 assert(D);
40
41 // Let U32(X) denote a uint32_t containing the upper 32 bits of a 64 bit
42 // integer X.
43 // Let L32(X) denote a uint32_t containing the lower 32 bits of a 64 bit
44 // integer X.
45 // Let X[A, B] with A <= B denote bits A through B of the integer X.
46 // Let (A | B) denote the concatination of two 32 bit ints, A and B.
47 // IOW X = (A | B) => U32(X) == A && L32(X) == B
48 //
49 // compute M = val * N (a 96 bit int)
50 // ---------------------------------
51 // tmp2 = U32(val) * N (a 64 bit int)
52 // tmp1 = L32(val) * N (a 64 bit int)
53 // which means
54 // M = val * N = (tmp2 << 32) + tmp1
55 tmp2 = (val >> 32) * N;
56 tmp1 = (val & UINT32_MAX) * N;
57
58 // compute M[32, 95]
59 // tmp2 = tmp2 + U32(tmp1)
60 // = (U32(val) * N) + U32(L32(val) * N)
61 // = M[32, 95]
62 tmp2 += tmp1 >> 32;
63
64 // if M[64, 95] >= D, then M/D has bits > 63 set and we have
65 // an overflow.
66 if ((tmp2 >> 32) >= D) {
67 *res = UINT64_MAX;
68 return false;
69 }
70
71 // Divide. Going in we know
72 // tmp2 = M[32, 95]
73 // U32(tmp2) < D
74 r = tmp2 % D;
75 tmp2 /= D;
76
77 // At this point
78 // tmp1 = L32(val) * N
79 // tmp2 = M[32, 95] / D
80 // = (M / D)[32, 95]
81 // r = M[32, 95] % D
82 // U32(tmp2) = 0
83 //
84 // compute tmp1 = (r | M[0, 31])
85 tmp1 = (tmp1 & UINT32_MAX) | ((uint64_t)r << 32);
86
87 // Divide again. Keep the remainder around in order to round properly.
88 r = tmp1 % D;
89 tmp1 /= D;
90
91 // At this point
92 // tmp2 = (M / D)[32, 95]
93 // tmp1 = (M / D)[ 0, 31]
94 // r = M % D
95 // U32(tmp1) = 0
96 // U32(tmp2) = 0
97
98 // Pack the result and deal with the round-up case (As well as the
99 // remote possiblility over overflow in such a case).
100 *res = (tmp2 << 32) | tmp1;
101 if (r && round_up_not_down) {
102 ++(*res);
103 if (!(*res)) {
104 *res = UINT64_MAX;
105 return false;
106 }
107 }
108
109 return true;
110 }
111
linear_transform_s64_to_s64(int64_t val,int64_t basis1,int32_t N,uint32_t D,bool invert_frac,int64_t basis2,int64_t * out)112 static bool linear_transform_s64_to_s64(
113 int64_t val,
114 int64_t basis1,
115 int32_t N,
116 uint32_t D,
117 bool invert_frac,
118 int64_t basis2,
119 int64_t* out) {
120 uint64_t scaled, res;
121 uint64_t abs_val;
122 bool is_neg;
123
124 if (!out)
125 return false;
126
127 // Compute abs(val - basis_64). Keep track of whether or not this delta
128 // will be negative after the scale opertaion.
129 if (val < basis1) {
130 is_neg = true;
131 abs_val = basis1 - val;
132 } else {
133 is_neg = false;
134 abs_val = val - basis1;
135 }
136
137 if (N < 0)
138 is_neg = !is_neg;
139
140 if (!scale_u64_to_u64(abs_val,
141 invert_frac ? D : ABS(N),
142 invert_frac ? ABS(N) : D,
143 &scaled,
144 is_neg))
145 return false; // overflow/undeflow
146
147 // if scaled is >= 0x8000<etc>, then we are going to overflow or
148 // underflow unless ABS(basis2) is large enough to pull us back into the
149 // non-overflow/underflow region.
150 if (scaled & INT64_MIN) {
151 if (is_neg && (basis2 < 0))
152 return false; // certain underflow
153
154 if (!is_neg && (basis2 >= 0))
155 return false; // certain overflow
156
157 if (ABS(basis2) <= static_cast<int64_t>(scaled & INT64_MAX))
158 return false; // not enough
159
160 // Looks like we are OK
161 *out = (is_neg ? (-scaled) : scaled) + basis2;
162 } else {
163 // Scaled fits within signed bounds, so we just need to check for
164 // over/underflow for two signed integers. Basically, if both scaled
165 // and basis2 have the same sign bit, and the result has a different
166 // sign bit, then we have under/overflow. An easy way to compute this
167 // is
168 // (scaled_signbit XNOR basis_signbit) &&
169 // (scaled_signbit XOR res_signbit)
170 // ==
171 // (scaled_signbit XOR basis_signbit XOR 1) &&
172 // (scaled_signbit XOR res_signbit)
173
174 if (is_neg)
175 scaled = -scaled;
176 res = scaled + basis2;
177
178 if ((scaled ^ basis2 ^ INT64_MIN) & (scaled ^ res) & INT64_MIN)
179 return false;
180
181 *out = res;
182 }
183
184 return true;
185 }
186
doForwardTransform(int64_t a_in,int64_t * b_out) const187 bool LinearTransform::doForwardTransform(int64_t a_in, int64_t* b_out) const {
188 if (0 == a_to_b_denom)
189 return false;
190
191 return linear_transform_s64_to_s64(a_in,
192 a_zero,
193 a_to_b_numer,
194 a_to_b_denom,
195 false,
196 b_zero,
197 b_out);
198 }
199
doReverseTransform(int64_t b_in,int64_t * a_out) const200 bool LinearTransform::doReverseTransform(int64_t b_in, int64_t* a_out) const {
201 if (0 == a_to_b_numer)
202 return false;
203
204 return linear_transform_s64_to_s64(b_in,
205 b_zero,
206 a_to_b_numer,
207 a_to_b_denom,
208 true,
209 a_zero,
210 a_out);
211 }
212
reduce(T * N,T * D)213 template <class T> void LinearTransform::reduce(T* N, T* D) {
214 T a, b;
215 if (!N || !D || !(*D)) {
216 assert(false);
217 return;
218 }
219
220 a = *N;
221 b = *D;
222
223 if (a == 0) {
224 *D = 1;
225 return;
226 }
227
228 // This implements Euclid's method to find GCD.
229 if (a < b) {
230 T tmp = a;
231 a = b;
232 b = tmp;
233 }
234
235 while (1) {
236 // a is now the greater of the two.
237 const T remainder = a % b;
238 if (remainder == 0) {
239 *N /= b;
240 *D /= b;
241 return;
242 }
243 // by swapping remainder and b, we are guaranteeing that a is
244 // still the greater of the two upon entrance to the loop.
245 a = b;
246 b = remainder;
247 }
248 };
249
250 template void LinearTransform::reduce<uint64_t>(uint64_t* N, uint64_t* D);
251 template void LinearTransform::reduce<uint32_t>(uint32_t* N, uint32_t* D);
252
reduce(int32_t * N,uint32_t * D)253 void LinearTransform::reduce(int32_t* N, uint32_t* D) {
254 if (N && D && *D) {
255 if (*N < 0) {
256 *N = -(*N);
257 reduce(reinterpret_cast<uint32_t*>(N), D);
258 *N = -(*N);
259 } else {
260 reduce(reinterpret_cast<uint32_t*>(N), D);
261 }
262 }
263 }
264
265 } // namespace android
266