fixed pt binary lib
I have the beginning of a fixed-pt binary lib. Fixed-pt is popular in DSP applications. It is integer arithmetic, but keeping track of the binary pt. It is built on top of constrained_value. It is small, so I include it here:
On Thu, 2008-07-31 at 07:47 -0400, Neal Becker wrote:
I have the beginning of a fixed-pt binary lib. Fixed-pt is popular in DSP applications. It is integer arithmetic, but keeping track of the binary pt.
It is built on top of constrained_value.
It is small, so I include it here:
I think you should return self& instead of self from most of the operators. Regards, Daniel
Daniel Frey wrote:
On Thu, 2008-07-31 at 07:47 -0400, Neal Becker wrote:
I have the beginning of a fixed-pt binary lib. Fixed-pt is popular in DSP applications. It is integer arithmetic, but keeping track of the binary pt.
It is built on top of constrained_value.
It is small, so I include it here:
I think you should return self& instead of self from most of the operators.
Yes, thanks (oversight)
Neal Becker wrote:
Daniel Frey wrote:
On Thu, 2008-07-31 at 07:47 -0400, Neal Becker wrote:
I have the beginning of a fixed-pt binary lib. Fixed-pt is popular in DSP applications. It is integer arithmetic, but keeping track of the binary pt.
It is built on top of constrained_value.
It is small, so I include it here: I think you should return self& instead of self from most of the operators.
Yes, thanks (oversight)
Be sure to check the size of your type with your compiler. Some embedded and DSP compilers in the past made types that derive from operators larger. Although I think that was partially due to a operators bug that was fixed a long time ago. There have been several fixed-point libraries discussed on this list you should search the archives including one very interesting one built on proto by Maurizio Vitale IIRC. Thanks, Michael Marcin
On Thu, 2008-07-31 at 12:16 -0500, Michael Marcin wrote:
Neal Becker wrote:
Daniel Frey wrote:
On Thu, 2008-07-31 at 07:47 -0400, Neal Becker wrote:
I have the beginning of a fixed-pt binary lib. Fixed-pt is popular in DSP applications. It is integer arithmetic, but keeping track of the binary pt.
It is built on top of constrained_value.
It is small, so I include it here: I think you should return self& instead of self from most of the operators.
Yes, thanks (oversight)
Be sure to check the size of your type with your compiler. Some embedded and DSP compilers in the past made types that derive from operators larger. Although I think that was partially due to a operators bug that was fixed a long time ago.
AFAIK that's not a bug in the operators library, but an inefficiency in the compiler's ABI (the EBO actually, which misses some optimization opportunities for multiply inheritance from empty base classes). If your compiler is inefficient, you can use the case class chaining technique offered/supported by the operators library (currently, Neal's code is *not* using it). See the documentation at <http://www.boost.org/doc/libs/1_35_0/libs/utility/operators.htm#chaining>. Regards, Daniel
Daniel Frey wrote:
On Thu, 2008-07-31 at 12:16 -0500, Michael Marcin wrote:
Be sure to check the size of your type with your compiler. Some embedded and DSP compilers in the past made types that derive from operators larger. Although I think that was partially due to a operators bug that was fixed a long time ago.
AFAIK that's not a bug in the operators library, but an inefficiency in the compiler's ABI (the EBO actually, which misses some optimization opportunities for multiply inheritance from empty base classes).
I was referring to http://svn.boost.org/trac/boost/ticket/979 which bit me on a project where we had a fixed-point 3 coordinate vector which had operators which contained a fixed-point types which had operators. Thanks, Michael Marcin
On Thu, 2008-07-31 at 14:20 -0500, Michael Marcin wrote:
Daniel Frey wrote:
On Thu, 2008-07-31 at 12:16 -0500, Michael Marcin wrote:
Be sure to check the size of your type with your compiler. Some embedded and DSP compilers in the past made types that derive from operators larger. Although I think that was partially due to a operators bug that was fixed a long time ago.
AFAIK that's not a bug in the operators library, but an inefficiency in the compiler's ABI (the EBO actually, which misses some optimization opportunities for multiply inheritance from empty base classes).
I was referring to http://svn.boost.org/trac/boost/ticket/979 which bit me on a project where we had a fixed-point 3 coordinate vector which had operators which contained a fixed-point types which had operators.
That won't help in Neal's code unless he starts using the base class chaining. Since empty_base only takes T into account (not U if available), Neal's code would still contain the same base class (empty< fixed_pt<int_bits,frac_bits,base_type,error_policy,round_policy> >) multiple times, which might lead to some bloat. Regards, Daniel
Daniel Frey wrote:
On Thu, 2008-07-31 at 14:20 -0500, Michael Marcin wrote:
Daniel Frey wrote:
On Thu, 2008-07-31 at 12:16 -0500, Michael Marcin wrote:
Be sure to check the size of your type with your compiler. Some embedded and DSP compilers in the past made types that derive from operators larger. Although I think that was partially due to a operators bug that was fixed a long time ago.
AFAIK that's not a bug in the operators library, but an inefficiency in the compiler's ABI (the EBO actually, which misses some optimization opportunities for multiply inheritance from empty base classes).
I was referring to http://svn.boost.org/trac/boost/ticket/979 which bit me on a project where we had a fixed-point 3 coordinate vector which had operators which contained a fixed-point types which had operators.
That won't help in Neal's code unless he starts using the base class chaining. Since empty_base only takes T into account (not U if available), Neal's code would still contain the same base class (empty< fixed_pt<int_bits,frac_bits,base_type,error_policy,round_policy> >) multiple times, which might lead to some bloat.
Regards, Daniel
I'll have to look into it. I'm using gcc-4.3.0, which I was guessing is not going to suffer bloat, so I didn't bother about using the chaining.
on Thu Jul 31 2008, Michael Marcin <mike.marcin-AT-gmail.com> wrote:
There have been several fixed-point libraries discussed on this list you should search the archives including one very interesting one built on proto by Maurizio Vitale IIRC.
Hey, yeah, cool idea. Then you can preserve full precision until the end of the full expression. -- Dave Abrahams BoostPro Computing http://www.boostpro.com
Neal Becker wrote:
I have the beginning of a fixed-pt binary lib. Fixed-pt is popular in DSP applications. It is integer arithmetic, but keeping track of the binary pt.
It is built on top of constrained_value.
It is small, so I include it here: Thanks, I am eagerly waiting for somebody to do something like this. Why did you choose to base it on constrained values? Couldn't I have used a bounded int as a base_type parameter, making your fixed_pt class not dependent on constrained_value anymore?
Also, do your implementations for /= and *= actually round? You actually seem to be using the round policy for conversions, is that right? Have you looked at the numeric/conversion stuff? (I haven't grokked it yet).
//self operator %=(base_type x) { val %= x; return *this; } What should this do? I think perhaps it should behave like fmod (?) --Johhn
John C. Femiani wrote:
Neal Becker wrote:
//self operator %=(base_type x) { val %= x; return *this; } What should this do? I think perhaps it should behave like fmod (?)
I had a fixed-point library with the purpose that code written should be able to choose to use fixed-point or floating point math at compile time without having to rewrite the algorithms. Instead of providing a % operator I provided an fmod function for parity with floating point types. /// Computes the remainder of fixed-point division /// @param x The numerator. /// @param y The divisor. /// @return the value x - i * y, for some integer i such that, if y is /// non-zero, the result has the same sign as x and magnitude /// less than the magnitude of y. template<std::size_t M, std::size_t F> inline fixed<M,F> fmod( fixed<M,F> x, fixed<M,F> y ) { return as_fixed(x.raw() % y.raw()); } /// shorthand for fmod( fixed, integer ) template<std::size_t M, std::size_t F, typename T > inline typename boost::enable_if<boost::is_integral<T>,fixed<M,F> >::type fmod( fixed<M,F> x, T y ) { return fmod( x, fixed<M,F>(y) ); } Thanks, Michael Marcin
"John C. Femiani" <john.femiani@asu.edu> writes:
Neal Becker wrote:
I have the beginning of a fixed-pt binary lib. Fixed-pt is popular in DSP applications. It is integer arithmetic, but keeping track of the binary pt.
Also, do your implementations for /= and *= actually round? You actually seem to be using the round policy for conversions, is that right? Have you looked at the numeric/conversion stuff? (I haven't grokked it yet).
The division throws away the lower "frac_bits" of the result: self operator/=(self const& x) { base_type tmp = (base_type)val / (base_type)x.val; val = tmp << (frac_bits); return *this; } I published a more fully-fledged fixed point type here a few months ago. You can get it from my website (along with a link to the DDJ article where I describe how it works): http://www.justsoftwaresolutions.co.uk/news/optimizing-applications-with-fix... It has a fixed number if integer and fractional bits (35.28). However, this allows it to support sin, cos, exp and log by using lookup tables. Anthony -- Anthony Williams | Just Software Solutions Ltd Custom Software Development | http://www.justsoftwaresolutions.co.uk Registered in England, Company Number 5478976. Registered Office: 15 Carrallack Mews, St Just, Cornwall, TR19 7UL
Anthony Williams wrote:
"John C. Femiani" <john.femiani@asu.edu> writes:
Also, do your implementations for /= and *= actually round? You actually seem to be using the round policy for conversions, is that right? Have you looked at the numeric/conversion stuff? (I haven't grokked it yet).
The division throws away the lower "frac_bits" of the result: So why doesn't it use the rounding policy? I kind of expect it to be used anytime quantization occurs.
--John
"John C. Femiani" <john.femiani@asu.edu> writes:
Anthony Williams wrote:
"John C. Femiani" <john.femiani@asu.edu> writes:
Also, do your implementations for /= and *= actually round? You actually seem to be using the round policy for conversions, is that right? Have you looked at the numeric/conversion stuff? (I haven't grokked it yet).
The division throws away the lower "frac_bits" of the result: So why doesn't it use the rounding policy? I kind of expect it to be used anytime quantization occurs.
I'll leave Neal to answer that (it's his code), but getting division right is a pain. If you followed the link to my code, you'll see that correct division is not trivial. Anthony -- Anthony Williams | Just Software Solutions Ltd Custom Software Development | http://www.justsoftwaresolutions.co.uk Registered in England, Company Number 5478976. Registered Office: 15 Carrallack Mews, St Just, Cornwall, TR19 7UL
Anthony Williams wrote:
"John C. Femiani" <john.femiani@asu.edu> writes:
Anthony Williams wrote:
"John C. Femiani" <john.femiani@asu.edu> writes:
Also, do your implementations for /= and *= actually round? You actually seem to be using the round policy for conversions, is that right? Have you looked at the numeric/conversion stuff? (I haven't grokked it yet).
The division throws away the lower "frac_bits" of the result: So why doesn't it use the rounding policy? I kind of expect it to be used anytime quantization occurs.
I'll leave Neal to answer that (it's his code), but getting division right is a pain. If you followed the link to my code, you'll see that correct division is not trivial.
Anthony
Multiplication is a bit easier to understand. I thought about mixed-arith mult, but I think it opens a can of worms. What do you want the result to be? So, I decided to only support arith of like-types. If you want the result of mult to be a larger type, then you simply ask for it: fixed<8,8> x, y; fixed<16,16> z = fixed<16,16> (x) * fixed<16,16> (y); As far as division, I wasn't quite sure and the implementation shown is only a guess. Suggestions welcome.
Neal Becker wrote:
Multiplication is a bit easier to understand. I thought about mixed-arith mult, but I think it opens a can of worms. What do you want the result to be? So, I decided to only support arith of like-types. If you want the result of mult to be a larger type, then you simply ask for it:
fixed<8,8> x, y; fixed<16,16> z = fixed<16,16> (x) * fixed<16,16> (y);
fixed<m0,f0> * fixed<m1,f1> naturally yields fixed<m0+m1,f0+f1> In my fixed point library I chose to make fixed<m,f>*fixed<m,f> return fixed<m,f> instead of fixed<2*m,2*f> but to make mixed type arithmetic return a detail type that is convertible to a fixed point type. namespace detail { template< typename PromotedType, std::size_t F0, std::size_t F1 > struct fixed_product { static const std::size_t frac_0 = F0; static const std::size_t frac_1 = F1; fixed_product( PromotedType r ) : data(r) {} PromotedType data; }; } // end namespace detail /// Product of 2 different fixed-point types /// @note this currently only supports the case we the two operands /// have the same bit resolution, it should be extended to /// cover the case of operands with different bit resolutions template< std::size_t M0, std::size_t F0, std::size_t M1, std::size_t F1 > inline typename boost::enable_if_c< M0+F0 == M1+F1, detail::fixed_product< typename detail::type_from_size<M0+F0>::next_size::value_type, F0, F1>
::type operator * ( fixed<M0,F0> a, fixed<M1,F1> b ) { typedef typename fixed<M0,F0>::next_type next_type; return detail::fixed_product<next_type,F0,F1>( static_cast<next_type>(a.raw()) * b.raw() ); }
template< std::size_t M, std::size_t F > template< std::size_t F0, std::size_t F1 > fixed<M,F>::fixed( detail::fixed_product<next_type,F0,F1> product , typename boost::enable_if< boost::mpl::greater_equal< boost::mpl::int_< F1+F0 > , boost::mpl::int_< F > > >::type* = 0 ) : data( static_cast<base_type>( product.data >> ((F1+F0)-F) ) ) { } template< std::size_t M, std::size_t F > template< std::size_t F0, std::size_t F1 > fixed<M,F>::fixed( detail::fixed_product<next_type,F0,F1> product , typename boost::enable_if< boost::mpl::less< boost::mpl::int_< F1+F0 > , boost::mpl::int_< F > > >::type* = 0 ) : data( static_cast<base_type>( product.data << (F-(F1+F0)) ) ) { } Thanks, Michael Marcin
Anthony Williams wrote:
"John C. Femiani" <john.femiani@asu.edu> writes:
Also, do your implementations for /= and *= actually round? You actually seem to be using the round policy for conversions, is that right? Have you looked at the numeric/conversion stuff? (I haven't grokked it yet).
The division throws away the lower "frac_bits" of the result:
Speaking of division, while I'd expect /= to return the same type as the LHS operand, the way I've previously implemented fixed point division is to return: result_int_bits = (numer_int_bits + denom_frac_bits + 1) result_frac_bits = (numer_frac_bits + denom_int_bits - 1) It's a little easier to understand the rationale, if you separately consider 1/x and x * y. For 1/x, returning the following avoids overflow and is reversible (if the result happens to be exact): result_int_bits = input_frac_bits + 1 result_frac_bits = denom_int_bits -1 The result of x * y is obviously: result_int_bits = x_int_bits + y_int_bits result_frac_bits = x_frac_bits + y_frac_bits Then, x/y is simply viewed as multiplication by the inverse of y, but computed in one shot. I haven't implemented /=, but I'd be tempted to compute it with less intermediate precision (i.e. only pre-scaling the numerator by 2^denom_frac_bits). I suppose, depending on rounding policy, you could compute it with an extra bit that gets rounded off. Matt
Matt Gruenke wrote:
Speaking of division, while I'd expect /= to return the same type as the LHS operand, the way I've previously implemented fixed point division is to return: result_int_bits = (numer_int_bits + denom_frac_bits + 1) result_frac_bits = (numer_frac_bits + denom_int_bits - 1) So then what happens to the value of X after I say Z = X /=Y? Do you mean '/', instead of '/='?
--John
John C. Femiani wrote:
Matt Gruenke wrote:
Speaking of division, while I'd expect /= to return the same type as the LHS operand, the way I've previously implemented fixed point division is to return: result_int_bits = (numer_int_bits + denom_frac_bits + 1) result_frac_bits = (numer_frac_bits + denom_int_bits - 1) So then what happens to the value of X after I say Z = X /=Y? Do you mean '/', instead of '/='?
Yes, I was describing the behavior of my operator/(). I hadn't written an operator/=(). Matt
Matt Gruenke wrote:
John C. Femiani wrote:
Matt Gruenke wrote:
Speaking of division, while I'd expect /= to return the same type as the LHS operand, the way I've previously implemented fixed point division is to return: result_int_bits = (numer_int_bits + denom_frac_bits + 1) result_frac_bits = (numer_frac_bits + denom_int_bits - 1) So then what happens to the value of X after I say Z = X /=Y? Do you mean '/', instead of '/='?
Yes, I was describing the behavior of my operator/(). I hadn't written an operator/=().
I think perhaps the issue with division is that rounding policy should be rethought. Also:
base_type tmp = (base_type)val / (base_type)x.val; val = tmp << (frac_bits); should be something like
static const int DENOMINATOR= 1 << frac_bits; base_type tmp = rounding_policy::divide(val * DENOMINATOR, x.val); I think the old way lost a lot of precision because it divided before multiplying.e.g. Then for multiplying:
base_type tmp = (base_type)val * (base_type)x.val; val = tmp >> (frac_bits);
could possibly be: val = rounding_policy::divide( val * x.val, DENOMINATOR); Another way to think of it is that the rounding policy should be a part of the base_type, decoupled from the fixed pt class. Then one can just assume that val * x.val / DENOMINATOR rounds the way the user intended. You can provide separate wrappers that change the bounds as well as the rounding behavior, and your code should be less coupled. To create a fixed point class, I could say: fixed_pt<8,24, bounded_int<int, MIN, MAX, true, true> > or fixed_pt<8,24, with_rounding< nearest, bounded_int<int, MIN, MAX, true, true> > > assuming that you had a template class called 'with_rounding' and a struct for the 'nearest' rounding policy. You can follow the trend of the constrained_value authors and make 'nearest' a function, so that ::boost::function can be used and the rounding policy can be changed at runtime. --John
John C. Femiani wrote:
Matt Gruenke wrote:
John C. Femiani wrote:
Matt Gruenke wrote:
Speaking of division, while I'd expect /= to return the same type as the LHS operand, the way I've previously implemented fixed point division is to return: result_int_bits = (numer_int_bits + denom_frac_bits + 1) result_frac_bits = (numer_frac_bits + denom_int_bits - 1) So then what happens to the value of X after I say Z = X /=Y? Do you mean '/', instead of '/='?
Yes, I was describing the behavior of my operator/(). I hadn't written an operator/=().
I think perhaps the issue with division is that rounding policy should be rethought. Also:
base_type tmp = (base_type)val / (base_type)x.val; val = tmp << (frac_bits); should be something like
static const int DENOMINATOR= 1 << frac_bits; base_type tmp = rounding_policy::divide(val * DENOMINATOR, x.val);
I think the old way lost a lot of precision because it divided before multiplying.e.g.
By replying to my message, it seems like you might be confusing two implementations. The one I was describing is my own private code, which worked the way you're suggesting.
Then for multiplying:
base_type tmp = (base_type)val * (base_type)x.val; val = tmp >> (frac_bits);
could possibly be:
val = rounding_policy::divide( val * x.val, DENOMINATOR);
The approach I took for operator*() was to return a larger type that lost no precision or range from the multiplication. Then the rounding happens when you assign the result to a smaller type. This both gives the user a natural way to control how much precision is preserved and avoids overhead from any unnecessary rounding (i.e. rounding more frequently than the user might want or need). Matt
Matt Gruenke wrote:
By replying to my message, it seems like you might be confusing two implementations. The one I was describing is my own private code, which worked the way you're suggesting.
Right, sorry for the ambiguity. I was referring to the code attached to the message by the OP. --John
On Fri, Aug 1, 2008 at 04:21, Matt Gruenke <mgruenke@intellivid.com> wrote:
The result of x * y is obviously: result_int_bits = x_int_bits + y_int_bits result_frac_bits = x_frac_bits + y_frac_bits
For unsigned, yes, but if int_bits includes the sign bit, it's result_int_bits = x_int_bits + y_int_bits - 1 (And it looked in the code like int_bits did, which does match what boost.integer does, and means you can do fixed_pt<16,16,int32_t>, which I think is more intuitive than that not actually fitting.)
participants (8)
-
Anthony Williams -
Daniel Frey -
David Abrahams -
John C. Femiani -
Matt Gruenke -
Michael Marcin -
Neal Becker -
Scott McMurray