Go to the documentation of this file.00001
00002 #include "vnl_rational.h"
00003
00004
00005
00006 #include <vnl/vnl_numeric_traits.h>
00007 #include <vcl_cassert.h>
00008
00009
00010
00011 vnl_rational::vnl_rational(double d)
00012 {
00013 bool sign = d<0;
00014 if (sign) d = -d;
00015
00016
00017 long den=0L, num=1L, prev_den=1L, prev_num=0L;
00018
00019 while (d*num < 1e9 && d*den < 1e9) {
00020 long a = (long)d;
00021 d -= a;
00022 long temp = num; num = a*num + prev_num; prev_num = temp;
00023 temp = den; den = a*den + prev_den; prev_den = temp;
00024 if (d < 1e-6) break;
00025 d = 1/d;
00026 }
00027 num_ = num; den_ = den;
00028 if (sign) num_ = -num_;
00029
00030 }
00031
00032
00033
00034
00035
00036 vnl_rational& vnl_rational::operator*=(vnl_rational const& r)
00037 {
00038 assert(num_!=0 || den_ != 0);
00039 long a = vnl_rational::gcd(r.numerator(),den_),
00040 b = vnl_rational::gcd(r.denominator(),num_);
00041 num_ /= b; den_ /= a;
00042 a = r.numerator()/a; b = r.denominator()/b;
00043
00044 double n = double(a) * double(num_),
00045 d = double(b) * double(den_);
00046 if (n < vnl_numeric_traits<long>::maxval && d < vnl_numeric_traits<long>::maxval)
00047 { num_ *= a; den_ *= b; normalize(); return *this; }
00048 else
00049 return *this = vnl_rational(n/d);
00050 }
00051
00052
00053
00054
00055 vnl_rational& vnl_rational::operator*=(long r)
00056 {
00057 long a = vnl_rational::gcd(r,den_);
00058 den_ /= a; r /= a;
00059
00060 double n = double(r) * double(num_);
00061 if (n < vnl_numeric_traits<long>::maxval)
00062 { num_ *= r; normalize(); return *this; }
00063 else
00064 return *this = vnl_rational(n/double(den_));
00065 }
00066
00067
00068
00069
00070
00071 vnl_rational& vnl_rational::operator/=(vnl_rational const& r)
00072 {
00073 assert(num_!=0 || den_ != 0);
00074 long a = vnl_rational::gcd(r.numerator(),num_),
00075 b = vnl_rational::gcd(r.denominator(),den_);
00076 num_ /= a; den_ /= b;
00077 a = r.numerator()/a; b = r.denominator()/b;
00078
00079 double n = double(b) * double(num_),
00080 d = double(a) * double(den_);
00081 if (n < vnl_numeric_traits<long>::maxval && d < vnl_numeric_traits<long>::maxval)
00082 { num_ *= b; den_ *= a; normalize(); return *this; }
00083 else
00084 return *this = vnl_rational(n/d);
00085 }
00086
00087
00088
00089
00090
00091 vnl_rational& vnl_rational::operator/=(long r)
00092 {
00093 assert(num_!=0 || r != 0);
00094 long a = vnl_rational::gcd(r,num_);
00095 num_ /= a; r /= a;
00096
00097 double d = double(r) * double(den_);
00098 if (d < vnl_numeric_traits<long>::maxval)
00099 { den_ *= r; normalize(); return *this; }
00100 else
00101 return *this = vnl_rational(double(num_)/d);
00102 }