Home | History | Annotate | Download | only in utils
      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 
     26 template<class T> static inline T ABS(T x) { return (x < 0) ? -x : x; }
     27 
     28 // Static math methods involving linear transformations
     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 
    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 
    187 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 
    200 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 
    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 
    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