libgfortran/io/write_float.def

   1 /* Copyright (C) 2007, 2008, 2009 Free Software Foundation, Inc.
   2    Contributed by Andy Vaught
   3    Write float code factoring to this file by Jerry DeLisle
   4    F2003 I/O support contributed by Jerry DeLisle
   5
   6 This file is part of the GNU Fortran 95 runtime library (libgfortran).
   7
   8 Libgfortran is free software; you can redistribute it and/or modify
   9 it under the terms of the GNU General Public License as published by
  10 the Free Software Foundation; either version 3, or (at your option)
  11 any later version.
  12
  13 Libgfortran is distributed in the hope that it will be useful,
  14 but WITHOUT ANY WARRANTY; without even the implied warranty of
  15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  16 GNU General Public License for more details.
  17
  18 Under Section 7 of GPL version 3, you are granted additional
  19 permissions described in the GCC Runtime Library Exception, version
  20 3.1, as published by the Free Software Foundation.
  21
  22 You should have received a copy of the GNU General Public License and
  23 a copy of the GCC Runtime Library Exception along with this program;
  24 see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
  25 <http://www.gnu.org/licenses/>.  */
  26
  27 #include "config.h"
  28
  29 typedef enum
  30 { S_NONE, S_MINUS, S_PLUS }
  31 sign_t;
  32
  33 /* Given a flag that indicates if a value is negative or not, return a
  34    sign_t that gives the sign that we need to produce.  */
  35
  36 static sign_t
  37 calculate_sign (st_parameter_dt *dtp, int negative_flag)
  38 {
  39   sign_t s = S_NONE;
  40
  41   if (negative_flag)
  42     s = S_MINUS;
  43   else
  44     switch (dtp->u.p.sign_status)
  45       {
  46       case SIGN_SP:     /* Show sign. */
  47         s = S_PLUS;
  48         break;
  49       case SIGN_SS:     /* Suppress sign. */
  50         s = S_NONE;
  51         break;
  52       case SIGN_S:      /* Processor defined. */
  53       case SIGN_UNSPECIFIED:
  54         s = options.optional_plus ? S_PLUS : S_NONE;
  55         break;
  56       }
  57
  58   return s;
  59 }
  60
  61
  62 /* Output a real number according to its format which is FMT_G free.  */
  63
  64 static void
  65 output_float (st_parameter_dt *dtp, const fnode *f, char *buffer, size_t size,
  66               int sign_bit, bool zero_flag, int ndigits, int edigits)
  67 {
  68   char *out;
  69   char *digits;
  70   int e;
  71   char expchar, rchar;
  72   format_token ft;
  73   int w;
  74   int d;
  75   /* Number of digits before the decimal point.  */
  76   int nbefore;
  77   /* Number of zeros after the decimal point.  */
  78   int nzero;
  79   /* Number of digits after the decimal point.  */
  80   int nafter;
  81   /* Number of zeros after the decimal point, whatever the precision.  */
  82   int nzero_real;
  83   int leadzero;
  84   int nblanks;
  85   int i;
  86   sign_t sign;
  87
  88   ft = f->format;
  89   w = f->u.real.w;
  90   d = f->u.real.d;
  91
  92   rchar = '5';
  93   nzero_real = -1;
  94
  95   /* We should always know the field width and precision.  */
  96   if (d < 0)
  97     internal_error (&dtp->common, "Unspecified precision");
  98
  99   sign = calculate_sign (dtp, sign_bit);
 100
 101   /* The following code checks the given string has punctuation in the correct
 102      places.  Uncomment if needed for debugging.
 103      if (d != 0 && ((buffer[2] != '.' && buffer[2] != ',')
 104                     || buffer[ndigits + 2] != 'e'))
 105        internal_error (&dtp->common, "printf is broken");  */
 106
 107   /* Read the exponent back in.  */
 108   e = atoi (&buffer[ndigits + 3]) + 1;
 109
 110   /* Make sure zero comes out as 0.0e0.   */
 111   if (zero_flag)
 112     {
 113       e = 0;
 114       if (compile_options.sign_zero == 1)
 115         sign = calculate_sign (dtp, sign_bit);
 116       else
 117         sign = calculate_sign (dtp, 0);
 118
 119       /* Handle special cases.  */
 120       if (w == 0)
 121         w = d + 2;
 122
 123       /* For this one we choose to not output a decimal point.
 124          F95 10.5.1.2.1  */
 125       if (w == 1 && ft == FMT_F)
 126         {
 127           out = write_block (dtp, w);
 128           if (out == NULL)
 129             return;
 130           *out = '0';
 131           return;
 132         }
 133
 134     }
 135
 136   /* Normalize the fractional component.  */
 137   buffer[2] = buffer[1];
 138   digits = &buffer[2];
 139
 140   /* Figure out where to place the decimal point.  */
 141   switch (ft)
 142     {
 143     case FMT_F:
 144       nbefore = e + dtp->u.p.scale_factor;
 145       if (nbefore < 0)
 146         {
 147           nzero = -nbefore;
 148           nzero_real = nzero;
 149           if (nzero > d)
 150             nzero = d;
 151           nafter = d - nzero;
 152           nbefore = 0;
 153         }
 154       else
 155         {
 156           nzero = 0;
 157           nafter = d;
 158         }
 159       expchar = 0;
 160       break;
 161
 162     case FMT_E:
 163     case FMT_D:
 164       i = dtp->u.p.scale_factor;
 165       if (d <= 0 && i == 0)
 166         {
 167           generate_error (&dtp->common, LIBERROR_FORMAT, "Precision not "
 168                           "greater than zero in format specifier 'E' or 'D'");
 169           return;
 170         }
 171       if (i <= -d || i >= d + 2)
 172         {
 173           generate_error (&dtp->common, LIBERROR_FORMAT, "Scale factor "
 174                           "out of range in format specifier 'E' or 'D'");
 175           return;
 176         }
 177
 178       if (!zero_flag)
 179         e -= i;
 180       if (i < 0)
 181         {
 182           nbefore = 0;
 183           nzero = -i;
 184           nafter = d + i;
 185         }
 186       else if (i > 0)
 187         {
 188           nbefore = i;
 189           nzero = 0;
 190           nafter = (d - i) + 1;
 191         }
 192       else /* i == 0 */
 193         {
 194           nbefore = 0;
 195           nzero = 0;
 196           nafter = d;
 197         }
 198
 199       if (ft == FMT_E)
 200         expchar = 'E';
 201       else
 202         expchar = 'D';
 203       break;
 204
 205     case FMT_EN:
 206       /* The exponent must be a multiple of three, with 1-3 digits before
 207          the decimal point.  */
 208       if (!zero_flag)
 209         e--;
 210       if (e >= 0)
 211         nbefore = e % 3;
 212       else
 213         {
 214           nbefore = (-e) % 3;
 215           if (nbefore != 0)
 216             nbefore = 3 - nbefore;
 217         }
 218       e -= nbefore;
 219       nbefore++;
 220       nzero = 0;
 221       nafter = d;
 222       expchar = 'E';
 223       break;
 224
 225     case FMT_ES:
 226       if (!zero_flag)
 227         e--;
 228       nbefore = 1;
 229       nzero = 0;
 230       nafter = d;
 231       expchar = 'E';
 232       break;
 233
 234     default:
 235       /* Should never happen.  */
 236       internal_error (&dtp->common, "Unexpected format token");
 237     }
 238
 239   /* Round the value.  The value being rounded is an unsigned magnitude.
 240      The ROUND_COMPATIBLE is rounding away from zero when there is a tie.  */
 241   switch (dtp->u.p.current_unit->round_status)
 242     {
 243       case ROUND_ZERO: /* Do nothing and truncation occurs.  */
 244         goto skip;
 245       case ROUND_UP:
 246         if (sign_bit)
 247           goto skip;
 248         rchar = '0';
 249         break;
 250       case ROUND_DOWN:
 251         if (!sign_bit)
 252           goto skip;
 253         rchar = '0';
 254         break;
 255       case ROUND_NEAREST:
 256         /* Round compatible unless there is a tie. A tie is a 5 with
 257            all trailing zero's.  */
 258         i = nafter + 1;
 259         if (digits[i] == '5')
 260           {
 261             for(i++ ; i < ndigits; i++)
 262               {
 263                 if (digits[i] != '0')
 264                   goto do_rnd;
 265               }
 266             /* It is a  tie so round to even.  */
 267             switch (digits[nafter])
 268               {
 269                 case '1':
 270                 case '3':
 271                 case '5':
 272                 case '7':
 273                 case '9':
 274                   /* If odd, round away from zero to even.  */
 275                   break;
 276                 default:
 277                   /* If even, skip rounding, truncate to even.  */
 278                   goto skip;
 279               }
 280           }
 281          /* Fall through.  */
 282       case ROUND_PROCDEFINED:
 283       case ROUND_UNSPECIFIED:
 284       case ROUND_COMPATIBLE:
 285         rchar = '5';
 286         /* Just fall through and do the actual rounding.  */
 287     }
 288
 289   do_rnd:
 290
 291   if (nbefore + nafter == 0)
 292     {
 293       ndigits = 0;
 294       if (nzero_real == d && digits[0] >= rchar)
 295         {
 296           /* We rounded to zero but shouldn't have */
 297           nzero--;
 298           nafter = 1;
 299           digits[0] = '1';
 300           ndigits = 1;
 301         }
 302     }
 303   else if (nbefore + nafter < ndigits)
 304     {
 305       ndigits = nbefore + nafter;
 306       i = ndigits;
 307       if (digits[i] >= rchar)
 308         {
 309           /* Propagate the carry.  */
 310           for (i--; i >= 0; i--)
 311             {
 312               if (digits[i] != '9')
 313                 {
 314                   digits[i]++;
 315                   break;
 316                 }
 317               digits[i] = '0';
 318             }
 319
 320           if (i < 0)
 321             {
 322               /* The carry overflowed.  Fortunately we have some spare
 323                  space at the start of the buffer.  We may discard some
 324                  digits, but this is ok because we already know they are
 325                  zero.  */
 326               digits--;
 327               digits[0] = '1';
 328               if (ft == FMT_F)
 329                 {
 330                   if (nzero > 0)
 331                     {
 332                       nzero--;
 333                       nafter++;
 334                     }
 335                   else
 336                     nbefore++;
 337                 }
 338               else if (ft == FMT_EN)
 339                 {
 340                   nbefore++;
 341                   if (nbefore == 4)
 342                     {
 343                       nbefore = 1;
 344                       e += 3;
 345                     }
 346                 }
 347               else
 348                 e++;
 349             }
 350         }
 351     }
 352
 353   skip:
 354
 355   /* Calculate the format of the exponent field.  */
 356   if (expchar)
 357     {
 358       edigits = 1;
 359       for (i = abs (e); i >= 10; i /= 10)
 360         edigits++;
 361
 362       if (f->u.real.e < 0)
 363         {
 364           /* Width not specified.  Must be no more than 3 digits.  */
 365           if (e > 999 || e < -999)
 366             edigits = -1;
 367           else
 368             {
 369               edigits = 4;
 370               if (e > 99 || e < -99)
 371                 expchar = ' ';
 372             }
 373         }
 374       else
 375         {
 376           /* Exponent width specified, check it is wide enough.  */
 377           if (edigits > f->u.real.e)
 378             edigits = -1;
 379           else
 380             edigits = f->u.real.e + 2;
 381         }
 382     }
 383   else
 384     edigits = 0;
 385
 386   /* Zero values always output as positive, even if the value was negative
 387      before rounding.  */
 388   for (i = 0; i < ndigits; i++)
 389     {
 390       if (digits[i] != '0')
 391         break;
 392     }
 393   if (i == ndigits)
 394     {
 395       /* The output is zero, so set the sign according to the sign bit unless
 396          -fno-sign-zero was specified.  */
 397       if (compile_options.sign_zero == 1)
 398         sign = calculate_sign (dtp, sign_bit);
 399       else
 400         sign = calculate_sign (dtp, 0);
 401     }
 402
 403   /* Pick a field size if none was specified.  */
 404   if (w <= 0)
 405     w = nbefore + nzero + nafter + (sign != S_NONE ? 2 : 1);
 406
 407   /* Work out how much padding is needed.  */
 408   nblanks = w - (nbefore + nzero + nafter + edigits + 1);
 409   if (sign != S_NONE)
 410     nblanks--;
 411
 412   if (dtp->u.p.g0_no_blanks)
 413     {
 414       w -= nblanks;
 415       nblanks = 0;
 416     }
 417
 418   /* Create the ouput buffer.  */
 419   out = write_block (dtp, w);
 420   if (out == NULL)
 421     return;
 422
 423   /* Check the value fits in the specified field width.  */
 424   if (nblanks < 0 || edigits == -1)
 425     {
 426       star_fill (out, w);
 427       return;
 428     }
 429
 430   /* See if we have space for a zero before the decimal point.  */
 431   if (nbefore == 0 && nblanks > 0)
 432     {
 433       leadzero = 1;
 434       nblanks--;
 435     }
 436   else
 437     leadzero = 0;
 438
 439   /* Pad to full field width.  */
 440
 441   if ( ( nblanks > 0 ) && !dtp->u.p.no_leading_blank)
 442     {
 443       memset (out, ' ', nblanks);
 444       out += nblanks;
 445     }
 446
 447   /* Output the initial sign (if any).  */
 448   if (sign == S_PLUS)
 449     *(out++) = '+';
 450   else if (sign == S_MINUS)
 451     *(out++) = '-';
 452
 453   /* Output an optional leading zero.  */
 454   if (leadzero)
 455     *(out++) = '0';
 456
 457   /* Output the part before the decimal point, padding with zeros.  */
 458   if (nbefore > 0)
 459     {
 460       if (nbefore > ndigits)
 461         {
 462           i = ndigits;
 463           memcpy (out, digits, i);
 464           ndigits = 0;
 465           while (i < nbefore)
 466             out[i++] = '0';
 467         }
 468       else
 469         {
 470           i = nbefore;
 471           memcpy (out, digits, i);
 472           ndigits -= i;
 473         }
 474
 475       digits += i;
 476       out += nbefore;
 477     }
 478
 479   /* Output the decimal point.  */
 480   *(out++) = dtp->u.p.current_unit->decimal_status == DECIMAL_POINT ? '.' : ',';
 481
 482   /* Output leading zeros after the decimal point.  */
 483   if (nzero > 0)
 484     {
 485       for (i = 0; i < nzero; i++)
 486         *(out++) = '0';
 487     }
 488
 489   /* Output digits after the decimal point, padding with zeros.  */
 490   if (nafter > 0)
 491     {
 492       if (nafter > ndigits)
 493         i = ndigits;
 494       else
 495         i = nafter;
 496
 497       memcpy (out, digits, i);
 498       while (i < nafter)
 499         out[i++] = '0';
 500
 501       digits += i;
 502       ndigits -= i;
 503       out += nafter;
 504     }
 505
 506   /* Output the exponent.  */
 507   if (expchar)
 508     {
 509       if (expchar != ' ')
 510         {
 511           *(out++) = expchar;
 512           edigits--;
 513         }
 514 #if HAVE_SNPRINTF
 515       snprintf (buffer, size, "%+0*d", edigits, e);
 516 #else
 517       sprintf (buffer, "%+0*d", edigits, e);
 518 #endif
 519       memcpy (out, buffer, edigits);
 520     }
 521
 522   if (dtp->u.p.no_leading_blank)
 523     {
 524       out += edigits;
 525       memset( out , ' ' , nblanks );
 526       dtp->u.p.no_leading_blank = 0;
 527     }
 528
 529 #undef STR
 530 #undef STR1
 531 #undef MIN_FIELD_WIDTH
 532 }
 533
 534
 535 /* Write "Infinite" or "Nan" as appropriate for the given format.  */
 536
 537 static void
 538 write_infnan (st_parameter_dt *dtp, const fnode *f, int isnan_flag, int sign_bit)
 539 {
 540   char * p, fin;
 541   int nb = 0;
 542
 543   if (f->format != FMT_B && f->format != FMT_O && f->format != FMT_Z)
 544     {
 545           nb =  f->u.real.w;
 546
 547           /* If the field width is zero, the processor must select a width
 548              not zero.  4 is chosen to allow output of '-Inf' or '+Inf' */
 549
 550           if (nb == 0) nb = 4;
 551           p = write_block (dtp, nb);
 552           if (p == NULL)
 553             return;
 554           if (nb < 3)
 555             {
 556               memset (p, '*',nb);
 557               return;
 558             }
 559
 560           memset(p, ' ', nb);
 561           if (!isnan_flag)
 562             {
 563               if (sign_bit)
 564                 {
 565
 566                   /* If the sign is negative and the width is 3, there is
 567                      insufficient room to output '-Inf', so output asterisks */
 568
 569                   if (nb == 3)
 570                     {
 571                       memset (p, '*',nb);
 572                       return;
 573                     }
 574
 575                   /* The negative sign is mandatory */
 576
 577                   fin = '-';
 578                 }
 579               else
 580
 581                   /* The positive sign is optional, but we output it for
 582                      consistency */
 583                   fin = '+';
 584
 585               if (nb > 8)
 586
 587                 /* We have room, so output 'Infinity' */
 588                 memcpy(p + nb - 8, "Infinity", 8);
 589               else
 590
 591                 /* For the case of width equals 8, there is not enough room
 592                    for the sign and 'Infinity' so we go with 'Inf' */
 593                 memcpy(p + nb - 3, "Inf", 3);
 594
 595               if (nb < 9 && nb > 3)
 596                 p[nb - 4] = fin;  /* Put the sign in front of Inf */
 597               else if (nb > 8)
 598                 p[nb - 9] = fin;  /* Put the sign in front of Infinity */
 599             }
 600           else
 601             memcpy(p + nb - 3, "NaN", 3);
 602           return;
 603         }
 604     }
 605
 606
 607 /* Returns the value of 10**d.  */
 608
 609 #define CALCULATE_EXP(x) \
 610 inline static GFC_REAL_ ## x \
 611 calculate_exp_ ## x  (int d)\
 612 {\
 613   int i;\
 614   GFC_REAL_ ## x r = 1.0;\
 615   for (i = 0; i< (d >= 0 ? d : -d); i++)\
 616     r *= 10;\
 617   r = (d >= 0) ? r : 1.0 / r;\
 618   return r;\
 619 }
 620
 621 CALCULATE_EXP(4)
 622
 623 CALCULATE_EXP(8)
 624
 625 #ifdef HAVE_GFC_REAL_10
 626 CALCULATE_EXP(10)
 627 #endif
 628
 629 #ifdef HAVE_GFC_REAL_16
 630 CALCULATE_EXP(16)
 631 #endif
 632 #undef CALCULATE_EXP
 633
 634 /* Generate corresponding I/O format for FMT_G and output.
 635    The rules to translate FMT_G to FMT_E or FMT_F from DEC fortran
 636    LRM (table 11-2, Chapter 11, "I/O Formatting", P11-25) is:
 637
 638    Data Magnitude                              Equivalent Conversion
 639    0< m < 0.1-0.5*10**(-d-1)                   Ew.d[Ee]
 640    m = 0                                       F(w-n).(d-1), n' '
 641    0.1-0.5*10**(-d-1)<= m < 1-0.5*10**(-d)     F(w-n).d, n' '
 642    1-0.5*10**(-d)<= m < 10-0.5*10**(-d+1)      F(w-n).(d-1), n' '
 643    10-0.5*10**(-d+1)<= m < 100-0.5*10**(-d+2)  F(w-n).(d-2), n' '
 644    ................                           ..........
 645    10**(d-1)-0.5*10**(-1)<= m <10**d-0.5       F(w-n).0,n(' ')
 646    m >= 10**d-0.5                              Ew.d[Ee]
 647
 648    notes: for Gw.d ,  n' ' means 4 blanks
 649           for Gw.dEe, n' ' means e+2 blanks  */
 650
 651 #define OUTPUT_FLOAT_FMT_G(x) \
 652 static void \
 653 output_float_FMT_G_ ## x (st_parameter_dt *dtp, const fnode *f, \
 654                       GFC_REAL_ ## x m, char *buffer, size_t size, \
 655                       int sign_bit, bool zero_flag, int ndigits, int edigits) \
 656 { \
 657   int e = f->u.real.e;\
 658   int d = f->u.real.d;\
 659   int w = f->u.real.w;\
 660   fnode *newf;\
 661   GFC_REAL_ ## x rexp_d;\
 662   int low, high, mid;\
 663   int ubound, lbound;\
 664   char *p;\
 665   int save_scale_factor, nb = 0;\
 666 \
 667   save_scale_factor = dtp->u.p.scale_factor;\
 668   newf = (fnode *) get_mem (sizeof (fnode));\
 669 \
 670   rexp_d = calculate_exp_ ## x (-d);\
 671   if ((m > 0.0 && m < 0.1 - 0.05 * rexp_d) || (rexp_d * (m + 0.5) >= 1.0) ||\
 672       ((m == 0.0) && !(compile_options.allow_std & GFC_STD_F2003)))\
 673     { \
 674       newf->format = FMT_E;\
 675       newf->u.real.w = w;\
 676       newf->u.real.d = d;\
 677       newf->u.real.e = e;\
 678       nb = 0;\
 679       goto finish;\
 680     }\
 681 \
 682   mid = 0;\
 683   low = 0;\
 684   high = d + 1;\
 685   lbound = 0;\
 686   ubound = d + 1;\
 687 \
 688   while (low <= high)\
 689     { \
 690       GFC_REAL_ ## x temp;\
 691       mid = (low + high) / 2;\
 692 \
 693       temp = (calculate_exp_ ## x (mid - 1) * (1 - 0.5 * rexp_d));\
 694 \
 695       if (m < temp)\
 696         { \
 697           ubound = mid;\
 698           if (ubound == lbound + 1)\
 699             break;\
 700           high = mid - 1;\
 701         }\
 702       else if (m > temp)\
 703         { \
 704           lbound = mid;\
 705           if (ubound == lbound + 1)\
 706             { \
 707               mid ++;\
 708               break;\
 709             }\
 710           low = mid + 1;\
 711         }\
 712       else\
 713         {\
 714           mid++;\
 715           break;\
 716         }\
 717     }\
 718 \
 719   if (e < 0)\
 720     nb = 4;\
 721   else\
 722     nb = e + 2;\
 723 \
 724   newf->format = FMT_F;\
 725   newf->u.real.w = f->u.real.w - nb;\
 726 \
 727   if (m == 0.0)\
 728     newf->u.real.d = d - 1;\
 729   else\
 730     newf->u.real.d = - (mid - d - 1);\
 731 \
 732   dtp->u.p.scale_factor = 0;\
 733 \
 734  finish:\
 735   output_float (dtp, newf, buffer, size, sign_bit, zero_flag, ndigits, \
 736                 edigits);\
 737   dtp->u.p.scale_factor = save_scale_factor;\
 738 \
 739   free_mem(newf);\
 740 \
 741   if (nb > 0 && !dtp->u.p.g0_no_blanks)\
 742     { \
 743       p = write_block (dtp, nb);\
 744       if (p == NULL)\
 745         return;\
 746       memset (p, ' ', nb);\
 747     }\
 748 }\
 749
 750 OUTPUT_FLOAT_FMT_G(4)
 751
 752 OUTPUT_FLOAT_FMT_G(8)
 753
 754 #ifdef HAVE_GFC_REAL_10
 755 OUTPUT_FLOAT_FMT_G(10)
 756 #endif
 757
 758 #ifdef HAVE_GFC_REAL_16
 759 OUTPUT_FLOAT_FMT_G(16)
 760 #endif
 761
 762 #undef OUTPUT_FLOAT_FMT_G
 763
 764
 765 /* Define a macro to build code for write_float.  */
 766
 767   /* Note: Before output_float is called, sprintf is used to print to buffer the
 768      number in the format +D.DDDDe+ddd. For an N digit exponent, this gives us
 769      (MIN_FIELD_WIDTH-5)-N digits after the decimal point, plus another one
 770      before the decimal point.
 771
 772      #   The result will always contain a decimal point, even if no
 773          digits follow it
 774
 775      -   The converted value is to be left adjusted on the field boundary
 776
 777      +   A sign (+ or -) always be placed before a number
 778
 779      MIN_FIELD_WIDTH  minimum field width
 780
 781      *   (ndigits-1) is used as the precision
 782
 783      e format: [-]d.dddeÂ±dd where there is one digit before the
 784        decimal-point character and the number of digits after it is
 785        equal to the precision. The exponent always contains at least two
 786        digits; if the value is zero, the exponent is 00.  */
 787
 788 #ifdef HAVE_SNPRINTF
 789
 790 #define DTOA \
 791 snprintf (buffer, size, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 792           "e", ndigits - 1, tmp);
 793
 794 #define DTOAL \
 795 snprintf (buffer, size, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 796           "Le", ndigits - 1, tmp);
 797
 798 #else
 799
 800 #define DTOA \
 801 sprintf (buffer, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 802          "e", ndigits - 1, tmp);
 803
 804 #define DTOAL \
 805 sprintf (buffer, "%+-#" STR(MIN_FIELD_WIDTH) ".*" \
 806          "Le", ndigits - 1, tmp);
 807
 808 #endif
 809
 810 #define WRITE_FLOAT(x,y)\
 811 {\
 812         GFC_REAL_ ## x tmp;\
 813         tmp = * (GFC_REAL_ ## x *)source;\
 814         sign_bit = signbit (tmp);\
 815         if (!isfinite (tmp))\
 816           { \
 817             write_infnan (dtp, f, isnan (tmp), sign_bit);\
 818             return;\
 819           }\
 820         tmp = sign_bit ? -tmp : tmp;\
 821         if (f->u.real.d == 0 && f->format == FMT_F\
 822             && dtp->u.p.scale_factor == 0)\
 823           {\
 824             if (tmp < 0.5)\
 825               tmp = 0.0;\
 826             else if (tmp < 1.0)\
 827               tmp = 1.0;\
 828           }\
 829         zero_flag = (tmp == 0.0);\
 830 \
 831         DTOA ## y\
 832 \
 833         if (f->format != FMT_G)\
 834           output_float (dtp, f, buffer, size, sign_bit, zero_flag, ndigits, \
 835                         edigits);\
 836         else \
 837           output_float_FMT_G_ ## x (dtp, f, tmp, buffer, size, sign_bit, \
 838                                     zero_flag, ndigits, edigits);\
 839 }\
 840
 841 /* Output a real number according to its format.  */
 842
 843 static void
 844 write_float (st_parameter_dt *dtp, const fnode *f, const char *source, int len)
 845 {
 846
 847 #if defined(HAVE_GFC_REAL_16) && __LDBL_DIG__ > 18
 848 # define MIN_FIELD_WIDTH 46
 849 #else
 850 # define MIN_FIELD_WIDTH 31
 851 #endif
 852 #define STR(x) STR1(x)
 853 #define STR1(x) #x
 854
 855   /* This must be large enough to accurately hold any value.  */
 856   char buffer[MIN_FIELD_WIDTH+1];
 857   int sign_bit, ndigits, edigits;
 858   bool zero_flag;
 859   size_t size;
 860
 861   size = MIN_FIELD_WIDTH+1;
 862
 863   /* printf pads blanks for us on the exponent so we just need it big enough
 864      to handle the largest number of exponent digits expected.  */
 865   edigits=4;
 866
 867   if (f->format == FMT_F || f->format == FMT_EN || f->format == FMT_G
 868       || ((f->format == FMT_D || f->format == FMT_E)
 869       && dtp->u.p.scale_factor != 0))
 870     {
 871       /* Always convert at full precision to avoid double rounding.  */
 872       ndigits = MIN_FIELD_WIDTH - 4 - edigits;
 873     }
 874   else
 875     {
 876       /* The number of digits is known, so let printf do the rounding.  */
 877       if (f->format == FMT_ES)
 878         ndigits = f->u.real.d + 1;
 879       else
 880         ndigits = f->u.real.d;
 881       if (ndigits > MIN_FIELD_WIDTH - 4 - edigits)
 882         ndigits = MIN_FIELD_WIDTH - 4 - edigits;
 883     }
 884
 885   switch (len)
 886     {
 887     case 4:
 888       WRITE_FLOAT(4,)
 889       break;
 890
 891     case 8:
 892       WRITE_FLOAT(8,)
 893       break;
 894
 895 #ifdef HAVE_GFC_REAL_10
 896     case 10:
 897       WRITE_FLOAT(10,L)
 898       break;
 899 #endif
 900 #ifdef HAVE_GFC_REAL_16
 901     case 16:
 902       WRITE_FLOAT(16,L)
 903       break;
 904 #endif
 905     default:
 906       internal_error (NULL, "bad real kind");
 907     }
 908 }