title: Rounding doubles to n decimal points using "round-to-even" in C++
body: <p>I want to implement a function in C++:</p>
  
  <pre><code>double round_to_even(double num, int decimal_places); 
  
  /* rounds the first argument to three decimal places 
  using round-to-even (unbiased rounding) */
  </code></pre>
  
  <p>which behaves the way rounding in R works, using a method known as round-to-even, unbiased rounding, or statistician's rounding: 
  <a href="https://stat.ethz.ch/R-manual/R-devel/library/base/html/Round.html" rel="nofollow">https://stat.ethz.ch/R-manual/R-devel/library/base/html/Round.html</a>
  The best way to illustrate this is by a few examples from R, where decimal_place=3:</p>
  
  <pre><code>&gt; round(0.1247,3) # rounding the first argument to three decimal places
  [1] 0.125
  
  &gt; round(0.1244,3)
  [1] 0.124
  </code></pre>
  
  <p>Notice how it behaves as you would expect, by rounding up when the digit to the right of the specified decimal_place is greater than 5 (here equal to 7), and truncating when it's smaller than 5 (here equal to 4). However, when the digit to the right of the specified decimal_place is equal to 5 (or at "midpoint") it behaves as follows: </p>
  
  <pre><code>&gt; round(0.1255,3) # will round up
  [1] 0.126
  
  &gt; round(.1275,3) # will round up
  [1] 0.128
  
  &gt; round(0.1265,3) # will truncate
  [1] 0.126
  
  &gt; round(0.1245,3) # will truncate
  [1] 0.124
  </code></pre>
  
  <p>by rounding up when the digit to the left of decimal_place is odd and truncating when even.</p>
  
  <p>Is there a concise way to implement this, and in particular without converting the digits to characters?</p>
  
  <p>Edit: 
  This is the best I could come up with:</p>
  
  <pre><code>double round_to_even(double number, int decimal_points) 
  /* rounds the first argument to three decimal places
  using round-to-even (unbiased rounding) */
  {   
      double num_left = number, num_right = number;
      int digit_left, digit_right;
  
      num_left *= pow(10, decimal_points + 1);
      digit_left = fmod(num_left, 10);
  
      if (digit_left == 5) 
      {   
          num_right *= pow(10, decimal_points);
          digit_right = fmod(num_right, 10);
  
          if (digit_right % 2 == 0) // if even 
              return floor(number * pow(10, decimal_points)) / pow(10, decimal_points);
          else // otherwise it's odd
              return ceil(number * pow(10, decimal_points)) / pow(10, decimal_points);
      }
      else { // standard round-to-nearest
          return round(number * pow(10, decimal_points)) / pow(10, decimal_points);
      }
  }
  </code></pre>
  
  <p>And I've tested it:</p>
  
  <pre><code>std::vector&lt;double&gt; test_vector({ 0.1247, 0.1244, 0.1255, 0.1275, 0.1265, 0.1245 });
  std::vector&lt;double&gt; expected_values({ 0.125, 0.124, 0.126, 0.128, 0.126, 0.124 });
  
  for (std::vector&lt;double&gt;::size_type i = 0; i &lt; test_vector.size(); i++)
      std::cout &lt;&lt; "expected: " &lt;&lt; expected_values[i] &lt;&lt; "\t got: " &lt;&lt; round_to_even(test_vector[i], 3) &lt;&lt; std::endl;
  </code></pre>
  
  <p>Which gives the following output:</p>
  
  <pre><code>expected: 0.125  got: 0.125 
  expected: 0.124  got: 0.124
  expected: 0.126  got: 0.126
  expected: 0.128  got: 0.128
  expected: 0.126  got: 0.126
  expected: 0.124  got: 0.124
  </code></pre>