src/if2img.cpp

   1 /*****************************************************************************
   2 **  This is part of the CTSim program
   3 **  Copyright (C) 1983-2000 Kevin Rosenberg
   4 **
   5 **  $Id: if2img.cpp,v 1.1 2000/06/07 07:55:30 kevin Exp $
   6 **  $Log: if2img.cpp,v $
   7 **  Revision 1.1  2000/06/07 07:55:30  kevin
   8 **  Initial import
   9 **
  10 **  Revision 1.1  2000/06/07 02:29:05  kevin
  11 **  Initial C++ versions
  12 **
  13 **  Revision 1.6  2000/06/03 07:57:51  kevin
  14 **  Fixed PNG 16-bit format
  15 **
  16 **  Revision 1.5  2000/06/03 06:29:57  kevin
  17 **  *** empty log message ***
  18 **
  19 **  Revision 1.4  2000/05/24 22:50:04  kevin
  20 **  Added support for new SGP library
  21 **
  22 **  Revision 1.3  2000/05/16 04:33:59  kevin
  23 **  Improved option processing
  24 **
  25 **  Revision 1.2  2000/05/08 20:02:32  kevin
  26 **  ANSI C changes
  27 **
  28 **  Revision 1.1.1.1  2000/04/28 13:02:44  kevin
  29 **  Initial CVS import for first public release
  30 **
  31 **
  32 **
  33 **  This program is free software; you can redistribute it and/or modify
  34 **  it under the terms of the GNU General Public License (version 2) as
  35 **  published by the Free Software Foundation.
  36 **
  37 **  This program is distributed in the hope that it will be useful,
  38 **  but WITHOUT ANY WARRANTY; without even the implied warranty of
  39 **  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  40 **  GNU General Public License for more details.
  41 **
  42 **  You should have received a copy of the GNU General Public License
  43 **  along with this program; if not, write to the Free Software
  44 **  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  45 ******************************************************************************/
  46
  47 /* FILE
  48  *   if2img.c           Convert an SDF file to a viewable format image
  49  */
  50
  51 #include "ct.h"
  52
  53 #if HAVE_PNG
  54 void sdf2d_to_png (ImageFile& im, char *outfile, int bitdepth, int nxcell, int nycell, double densmin, double densmax);
  55 #endif
  56 #if HAVE_GIF
  57 void sdf2d_to_gif (ImageFile& im, char *outfile, int nxcell, int nycell, double densmin, double densmax);
  58 #endif
  59 void sdf2d_to_pgm (ImageFile& im, char *outfile, int nxcell, int nycell, double densmin, double densmax);
  60 void sdf2d_to_pgmasc (ImageFile& im, char *outfile, int nxcell, int nycell, double densmin, double densmax);
  61
  62 enum { O_SCALE, O_MIN, O_MAX, O_AUTO, O_CENTER, O_STATS, O_FORMAT, O_LABELS,
  63        O_HELP, O_VERBOSE, O_VERSION, O_DEBUG };
  64
  65 static struct option my_options[] =
  66 {
  67   {"scale", 1, 0, O_SCALE},
  68   {"min", 1, 0, O_MIN},
  69   {"max", 1, 0, O_MAX},
  70   {"auto", 1, 0, O_AUTO},
  71   {"center", 1, 0, O_CENTER},
  72   {"format", 1, 0, O_FORMAT},
  73   {"labels", 0, 0, O_LABELS},
  74   {"stats", 0, 0, O_STATS},
  75   {"verbose", 0, 0, O_VERBOSE},
  76   {"debug", 0, 0, O_DEBUG},
  77   {"help", 0, 0, O_HELP},
  78   {"version", 0, 0, O_VERSION},
  79   {0, 0, 0, 0}
  80 };
  81
  82 enum { O_AUTO_FULL, O_AUTO_STD0_1, O_AUTO_STD0_5, O_AUTO_STD1, O_AUTO_STD2, O_AUTO_STD3 };
  83 #define O_AUTO_FULL_STR  "full"
  84 #define O_AUTO_STD0_1_STR  "std0.1"
  85 #define O_AUTO_STD0_5_STR  "std0.5"
  86 #define O_AUTO_STD1_STR  "std1"
  87 #define O_AUTO_STD2_STR  "std2"
  88 #define O_AUTO_STD3_STR  "std3"
  89
  90 enum { O_CENTER_MEAN, O_CENTER_MODE };
  91 #define O_CENTER_MEAN_STR  "mean"
  92 #define O_CENTER_MODE_STR  "mode"
  93
  94 enum { O_FORMAT_GIF, O_FORMAT_PNG, O_FORMAT_PNG16, O_FORMAT_PGM, O_FORMAT_PGMASC, O_FORMAT_DISP };
  95 #define O_FORMAT_GIF_STR   "gif"
  96 #define O_FORMAT_PNG_STR   "png"
  97 #define O_FORMAT_PNG16_STR   "png16"
  98 #define O_FORMAT_PGM_STR   "pgm"
  99 #define O_FORMAT_PGMASC_STR "pgmasc"
 100 #define O_FORMAT_DISP_STR    "disp"
 101
 102 void
 103 if2img_usage (const char *program)
 104 {
 105   fprintf(stdout, "usage: %s sdfname outfile [OPTIONS]\n", kbasename(program));
 106   fprintf(stdout, "Convert IF file to an image file\n");
 107   fprintf(stdout, "\n");
 108   fprintf(stdout, "     sdfname    Name of input SDF file\n");
 109   fprintf(stdout, "     outfile    Name of output file\n");
 110   fprintf(stdout, "     --format   Output format\n");
 111   fprintf(stdout, "         pgm    PGM (portable graymap) format (default)\n");
 112   fprintf(stdout, "         pgmasc PGM (portable graymap) ASCII format\n");
 113 #ifdef HAVE_PNG
 114   fprintf(stdout, "         png    PNG (8-bit) format\n");
 115   fprintf(stdout, "         png16  PNG (16-bit) format\n");
 116 #endif
 117 #ifdef HAVE_GIF
 118   fprintf(stdout, "         gif    GIF format\n");
 119 #endif
 120   fprintf(stdout, "         disp   Display on screen\n");
 121   fprintf(stdout, "     --center   Center of window\n");
 122   fprintf(stdout, "         mode   Mode is center of window (default)\n");
 123   fprintf(stdout, "         mean   Mean is center of window\n");
 124   fprintf(stdout, "     --auto     Set auto window\n");
 125   fprintf(stdout, "         full   Use full window (default)\n");
 126   fprintf(stdout, "         std0.1 Use 0.1 standard deviation about center\n");
 127   fprintf(stdout, "         std0.5 Use 0.5 standard deviation about center\n");
 128   fprintf(stdout, "         std1   Use one standard deviation about center\n");
 129   fprintf(stdout, "         std2   Use two standard deviations about center\n");
 130   fprintf(stdout, "         std3   Use three standard deviations about center\n");
 131   fprintf(stdout, "     --scale    Scaling factor for output size\n");
 132   fprintf(stdout, "     --min      Set minimum intensity\n");
 133   fprintf(stdout, "     --max      Set maximum intensity\n");
 134   fprintf(stdout, "     --stats    Print image statistics\n");
 135   fprintf(stdout, "     --labels   Print image labels\n");
 136   fprintf(stdout, "     --debug    Set debug mode\n");
 137   fprintf(stdout, "     --verbose  Set verbose mode\n");
 138   fprintf(stdout, "     --version  Print version\n");
 139   fprintf(stdout, "     --help     Print this help message\n");
 140 }
 141
 142
 143 int
 144 if2img_main (int argc, char *const argv[])
 145 {
 146   ImageFile* pim = NULL;
 147
 148   double densmin = HUGE_VAL, densmax = -HUGE_VAL;
 149   char *in_file, *out_file;
 150   int opt_verbose = 0;
 151   int opt_scale = 1;
 152   int opt_auto = O_AUTO_FULL;
 153   int opt_set_max = 0;
 154   int opt_set_min = 0;
 155   int opt_stats = 0;
 156   int opt_debug = 0;
 157   int opt_center = O_CENTER_MODE;
 158   int opt_format = O_FORMAT_PGM;
 159   int opt_labels = 0;
 160
 161   while (1)
 162     {
 163       int c = getopt_long (argc, argv, "", my_options, NULL);
 164       char *endptr = NULL;
 165       char *endstr;
 166
 167       if (c == -1)
 168         break;
 169
 170       switch (c)
 171         {
 172         case O_MIN:
 173           opt_set_min = 1;
 174           densmin = strtod(optarg, &endptr);
 175           endstr = optarg + strlen(optarg);
 176           if (endptr != endstr)
 177             {
 178               fprintf(stderr, "Error setting --min to %s\n", optarg);
 179               if2img_usage(argv[0]);
 180               return (1);
 181             }
 182           break;
 183         case O_MAX:
 184           opt_set_max = 1;
 185           densmax = strtod(optarg, &endptr);
 186           endstr = optarg + strlen(optarg);
 187           if (endptr != endstr)
 188             {
 189               fprintf(stderr, "Error setting --max to %s\n", optarg);
 190               if2img_usage(argv[0]);
 191               return (1);
 192             }
 193           break;
 194         case O_SCALE:
 195           opt_scale = strtol(optarg, &endptr, 10);
 196           endstr = optarg + strlen(optarg);
 197           if (endptr != endstr)
 198             {
 199               fprintf(stderr,"Error setting --scale to %s\n", optarg);
 200               if2img_usage(argv[0]);
 201               return (1);
 202             }
 203           break;
 204         case O_AUTO:
 205           if (strcmp(optarg, O_AUTO_FULL_STR) == 0)
 206             opt_auto = O_AUTO_FULL;
 207           else if (strcmp(optarg, O_AUTO_STD1_STR) == 0)
 208             opt_auto = O_AUTO_STD1;
 209           else if (strcmp(optarg, O_AUTO_STD0_5_STR) == 0)
 210             opt_auto = O_AUTO_STD0_5;
 211           else if (strcmp(optarg, O_AUTO_STD0_1_STR) == 0)
 212             opt_auto = O_AUTO_STD0_1;
 213           else if (strcmp(optarg, O_AUTO_STD2_STR) == 0)
 214             opt_auto = O_AUTO_STD2;
 215           else if (strcmp(optarg, O_AUTO_STD3_STR) == 0)
 216             opt_auto = O_AUTO_STD3;
 217           else
 218             {
 219               fprintf(stderr, "Invalid auto mode %s\n", optarg);
 220               if2img_usage(argv[0]);
 221               return (1);
 222             }
 223                 break;
 224         case O_CENTER:
 225           if (strcmp(optarg, O_CENTER_MEAN_STR) == 0)
 226             opt_center = O_CENTER_MEAN;
 227           else if (strcmp(optarg, O_CENTER_MODE_STR) == 0)
 228             opt_center = O_CENTER_MODE;
 229           else
 230             {
 231               fprintf(stderr, "Invalid center mode %s\n", optarg);
 232               if2img_usage(argv[0]);
 233               return (1);
 234             }
 235                 break;
 236         case O_FORMAT:
 237           if (strcmp(optarg, O_FORMAT_PGM_STR) == 0)
 238             opt_format = O_FORMAT_PGM;
 239           else if (strcmp(optarg, O_FORMAT_PGMASC_STR) == 0)
 240             opt_format = O_FORMAT_PGMASC;
 241 #if HAVE_PNG
 242           else if (strcmp(optarg, O_FORMAT_PNG_STR) == 0)
 243             opt_format = O_FORMAT_PNG;
 244           else if (strcmp(optarg, O_FORMAT_PNG16_STR) == 0)
 245             opt_format = O_FORMAT_PNG16;
 246 #endif
 247 #if HAVE_GIF
 248           else if (strcmp(optarg, O_FORMAT_GIF_STR) == 0)
 249             opt_format = O_FORMAT_GIF;
 250 #endif
 251           else if (strcmp(optarg, O_FORMAT_DISP_STR) == 0)
 252             opt_format = O_FORMAT_DISP;
 253           else {
 254               fprintf(stderr, "Invalid format mode %s\n", optarg);
 255               if2img_usage(argv[0]);
 256               return (1);
 257             }
 258           break;
 259         case O_LABELS:
 260           opt_labels = 1;
 261           break;
 262         case O_VERBOSE:
 263           opt_verbose = 1;
 264           break;
 265         case O_DEBUG:
 266           opt_debug = 1;
 267           break;
 268         case O_STATS:
 269           opt_stats = 1;
 270           break;
 271         case O_VERSION:
 272 #ifdef VERSION
 273           fprintf(stdout, "Version %s\n", VERSION);
 274 #else
 275           fprintf(stderr, "Unknown version number");
 276 #endif
 277           exit(0);
 278         case O_HELP:
 279         case '?':
 280           if2img_usage(argv[0]);
 281           return (0);
 282         default:
 283           if2img_usage(argv[0]);
 284           return (1);
 285         }
 286     }
 287
 288   if ((opt_format == O_FORMAT_DISP && optind + 1 != argc)
 289       || (opt_format != O_FORMAT_DISP && optind + 2 != argc))
 290     {
 291         if2img_usage(argv[0]);
 292         return (1);
 293     }
 294
 295   in_file = argv[optind];
 296
 297   if (opt_format != O_FORMAT_DISP)
 298       out_file = argv[optind+1];
 299   else out_file = NULL;
 300
 301   pim = new ImageFile (in_file);
 302   ImageFile& im = *pim;
 303   if (! im.adf.fileRead()) {
 304     fprintf(stderr, "File %s does not exist\n", in_file);
 305     return (1);
 306   }
 307
 308   if (opt_labels) {
 309     int nlabels = im.adf.getNumLabels();
 310
 311     for (int i = 0; i < nlabels; i++) {
 312       Array2dFileLabel label;
 313       im.adf.labelRead (label, i);
 314       string str;
 315       label.getDateString (str);
 316
 317       if (label.getLabelType() == Array2dFileLabel::L_HISTORY) {
 318         cout << "History: " << label.getLabelString() << endl;
 319         cout << "  calc time = " << label.getCalcTime() << " secs" << endl;
 320         cout << "  Timestamp = " << str << endl;
 321       } else if (label.getLabelType() == Array2dFileLabel::L_USER) {
 322         cout << "Note: " <<  label.getLabelString() << endl;
 323         cout << "  Timestamp = %s" << str << endl;
 324       }
 325     }
 326   }
 327
 328   if (opt_stats || (! (opt_set_max && opt_set_min))) {
 329     double minfound = HUGE_VAL, maxfound = -HUGE_VAL;
 330     double mode = 0, mean = 0, stddev = 0, window = 0;
 331     double spread;
 332     int hist[256];
 333     int ibin, nbin = 256;
 334     int max_bin, max_binindex;
 335     double maxbin;
 336     int nx = im.nx();
 337     int ny = im.ny();
 338     ImageFileArray v = im.getArray();
 339
 340     maxbin = nbin - 1;
 341     for (int ix = 0; ix < nx; ix++) {
 342       for (int iy = 0; iy < ny; iy++) {
 343         if (v[ix][iy] > maxfound)
 344           maxfound = v[ix][iy];
 345         if (v[ix][iy] < minfound)
 346           minfound = v[ix][iy];
 347         mean += v[ix][iy];
 348       }
 349     }
 350     spread = maxfound - minfound;
 351     if (spread == 0)
 352       mode = minfound;
 353     else {
 354       for (ibin = 0; ibin < nbin; ibin++)
 355         hist[ibin] = 0;
 356       for (int ix = 0; ix < nx; ix++) {
 357         for (int iy = 0; iy < ny; iy++) {
 358           int b = (int) ((((v[ix][iy] - minfound) / spread) * (double) maxbin) + 0.5);
 359           hist[b]++;
 360         }
 361       }
 362       max_binindex = 0;
 363       max_bin = -1;
 364       for (ibin = 0; ibin < nbin; ibin++) {
 365         if (hist[ibin] > max_bin) {
 366           max_bin = hist[ibin];
 367           max_binindex = ibin;
 368         }
 369       }
 370       mode = (((double) max_binindex) * spread / ((double) maxbin)) + minfound;
 371     }
 372
 373     if (opt_auto == O_AUTO_FULL) {
 374         if (! opt_set_max)
 375           densmax = maxfound;
 376         if (! opt_set_min)
 377           densmin = minfound;
 378     }
 379     if (opt_stats || opt_auto != O_AUTO_FULL) {
 380       mean /= nx * ny;
 381       for (int ix = 0; ix < nx; ix++) {
 382         for (int iy = 0; iy < ny; iy++) {
 383             double diff = (v[ix][iy] - mean);
 384             stddev += diff * diff;
 385         }
 386       }
 387       stddev = sqrt(stddev / (nx * ny));
 388       if (opt_auto == O_AUTO_FULL)
 389         ;
 390       else if (opt_auto == O_AUTO_STD1)
 391             window = stddev;
 392       else if (opt_auto == O_AUTO_STD0_1)
 393         window = stddev * 0.1;
 394       else if (opt_auto == O_AUTO_STD0_5)
 395         window = stddev * 0.5;
 396           else if (opt_auto == O_AUTO_STD2)
 397             window = stddev * 2;
 398       else if (opt_auto == O_AUTO_STD3)
 399         window = stddev * 3;
 400       else {
 401         fprintf(stderr, "Internal Error: Invalid auto mode %d\n", opt_auto);
 402         return (1);
 403       }
 404     }
 405     if (opt_stats) {
 406       fprintf(stdout,"nx=%d\n", nx);
 407       fprintf(stdout,"ny=%d\n", ny);
 408       fprintf(stdout,"min=%f\n", minfound);
 409       fprintf(stdout,"max=%f\n", maxfound);
 410       fprintf(stdout,"mean=%f\n", mean);
 411       fprintf(stdout,"mode=%f\n", mode);
 412       fprintf(stdout,"stddev=%f\n", stddev);
 413     }
 414     if (opt_auto != O_AUTO_FULL) {
 415       double center;
 416
 417       if (opt_center == O_CENTER_MODE)
 418         center = mode;
 419       else if (opt_center == O_CENTER_MEAN)
 420         center = mean;
 421       else {
 422         fprintf(stderr, "Internal Error: Invalid center mode %d\n", opt_center);
 423         return (1);
 424       }
 425       if (! opt_set_max)
 426         densmax = center + window;
 427       if (! opt_set_min)
 428         densmin = center - window;
 429       }
 430   }
 431
 432   if (opt_stats) {
 433     fprintf(stdout,"min_disp=%f\n", densmin);
 434     fprintf(stdout,"max_disp=%f\n", densmax);
 435   }
 436
 437   if (opt_format == O_FORMAT_PGM)
 438     sdf2d_to_pgm (im, out_file, opt_scale, opt_scale, densmin, densmax);
 439   else if (opt_format == O_FORMAT_PGMASC)
 440     sdf2d_to_pgmasc (im, out_file, opt_scale, opt_scale, densmin, densmax);
 441 #if HAVE_PNG
 442   else if (opt_format == O_FORMAT_PNG)
 443     sdf2d_to_png (im, out_file, 8, opt_scale, opt_scale, densmin, densmax);
 444   else if (opt_format == O_FORMAT_PNG16)
 445     sdf2d_to_png (im, out_file, 16, opt_scale, opt_scale, densmin, densmax);
 446 #endif
 447 #if HAVE_GIF
 448   else if (opt_format == O_FORMAT_GIF)
 449     sdf2d_to_gif (im, out_file, opt_scale, opt_scale, densmin, densmax);
 450 #endif
 451   else if (opt_format == O_FORMAT_DISP) {
 452 #if HAVE_SGP
 453     // image_display_scale (im, opt_scale, densmin, densmax);
 454     //  cio_kb_getc();
 455       sgp2_close(sgp2_get_active_win());
 456 #endif
 457   }
 458   else
 459     {
 460       fprintf(stderr, "Internal Error: Invalid format mode %d\n", opt_format);
 461       return (1);
 462     }
 463   return (0);
 464 }
 465
 466
 467 void
 468 sdf2d_to_pgm (ImageFile& im, char *outfile, int nxcell, int nycell, double densmin, double densmax)
 469 {
 470   FILE *fp;
 471   int irow;
 472   unsigned char *rowp;
 473   int nx = im.nx();
 474   int ny = im.ny();
 475   ImageFileArray v = im.getArray();
 476
 477   rowp = (unsigned char *) sys_alloc(nx * nxcell,"sdf2d_to_img");
 478   if (rowp == NULL)
 479     return;
 480
 481   if ((fp = fopen (outfile, "wb")) == NULL)
 482      return;
 483
 484   fprintf(fp, "P5\n");
 485   fprintf(fp, "%d %d\n", nx, ny);
 486   fprintf(fp, "255\n");
 487
 488   for (irow = ny - 1; irow >= 0; irow--)
 489     {
 490       int icol, ir;
 491
 492       for (icol = 0; icol < nx; icol++)
 493         {
 494           double dens;
 495           int p;
 496           int pos = icol * nxcell;
 497           dens = (v[icol][irow] - densmin) / (densmax - densmin);
 498           if (dens < 0)
 499             dens = 0;
 500           else if (dens > 1)
 501             dens = 1;
 502           for (p = pos; p < pos + nxcell; p++)
 503             {
 504               rowp[p] = (unsigned int) (dens * 255.);
 505             }
 506         }
 507       for (ir = 0; ir < nycell; ir++)
 508         {
 509           int ic;
 510           for (ic = 0; ic < nx * nxcell; ic++)
 511             fprintf(fp, "%c ", rowp[ic]);
 512         }
 513     }
 514   sys_free(rowp, "if2img");
 515
 516   fclose(fp);
 517 }
 518
 519 void
 520 sdf2d_to_pgmasc (ImageFile& im, char *outfile, int nxcell, int nycell, double densmin, double densmax)
 521 {
 522   FILE *fp;
 523   int irow;
 524   unsigned char *rowp;
 525   int nx = im.nx();
 526   int ny = im.ny();
 527   ImageFileArray v = im.getArray();
 528
 529   rowp = (unsigned char *) sys_alloc(nx * nxcell,"sdf2d_to_img");
 530   if (rowp == NULL)
 531     return;
 532
 533   if ((fp = fopen (outfile, "wb")) == NULL)
 534      return;
 535
 536   fprintf(fp, "P2\n");
 537   fprintf(fp, "%d %d\n", nx, ny);
 538   fprintf(fp, "255\n");
 539
 540   for (irow = ny - 1; irow >= 0; irow--)
 541     {
 542       int icol, ir;
 543
 544       for (icol = 0; icol < nx; icol++)
 545         {
 546           double dens;
 547           int p;
 548           int pos = icol * nxcell;
 549           dens = (v[icol][irow] - densmin) / (densmax - densmin);
 550           if (dens < 0)
 551             dens = 0;
 552           else if (dens > 1)
 553             dens = 1;
 554           for (p = pos; p < pos + nxcell; p++)
 555             {
 556               rowp[p] = (unsigned int) (dens * 255.);
 557             }
 558         }
 559       for (ir = 0; ir < nycell; ir++)
 560         {
 561           int ic;
 562           for (ic = 0; ic < nx * nxcell; ic++)
 563             fprintf(fp, "%d ", rowp[ic]);
 564           fprintf(fp, "\n");
 565         }
 566     }
 567   sys_free(rowp, "if2img");
 568
 569   fclose(fp);
 570 }
 571
 572
 573 #ifdef HAVE_PNG
 574 void
 575 sdf2d_to_png (ImageFile& im, char *outfile, int bitdepth, int nxcell, int nycell, double densmin, double densmax)
 576 {
 577   FILE *fp;
 578   png_structp png_ptr;
 579   png_infop info_ptr;
 580   int irow;
 581   unsigned char *rowp;
 582   double max_out_level = (1 << bitdepth) - 1;
 583   int nx = im.nx();
 584   int ny = im.ny();
 585   ImageFileArray v = im.getArray();
 586
 587   rowp = (unsigned char *) sys_alloc(nx * nxcell * (bitdepth / 8),"sdf2d_to_img");
 588   if (rowp == NULL)
 589     return;
 590
 591   if ((fp = fopen (outfile, "wb")) == NULL)
 592      return;
 593
 594   png_ptr = png_create_write_struct (PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);
 595   if (! png_ptr)
 596     return;
 597
 598   info_ptr = png_create_info_struct(png_ptr);
 599   if (! info_ptr)
 600     {
 601       png_destroy_write_struct(&png_ptr, (png_infopp) NULL);
 602       fclose(fp);
 603       return;
 604     }
 605
 606   if (setjmp(png_ptr->jmpbuf))
 607     {
 608       png_destroy_write_struct(&png_ptr, &info_ptr);
 609       fclose(fp);
 610       return;
 611     }
 612
 613   png_init_io(png_ptr, fp);
 614
 615   png_set_IHDR(png_ptr, info_ptr, nx * nxcell, ny * nycell, bitdepth, PNG_COLOR_TYPE_GRAY, PNG_INTERLACE_NONE, PNG_COMPRESSION_TYPE_BASE, PNG_FILTER_TYPE_DEFAULT);
 616
 617   png_write_info(png_ptr, info_ptr);
 618   for (irow = ny - 1; irow >= 0; irow--)
 619     {
 620       png_bytep row_pointer = rowp;
 621       int icol, ir;
 622
 623       for (icol = 0; icol < nx; icol++)
 624         {
 625           double dens;
 626           int p;
 627           unsigned int outval;
 628           int pos = icol * nxcell;
 629           dens = (v[icol][irow] - densmin) / (densmax - densmin);
 630           if (dens < 0)
 631             dens = 0;
 632           else if (dens > 1)
 633             dens = 1;
 634           outval = (unsigned int) (dens * max_out_level);
 635           for (p = pos; p < pos + nxcell; p++)
 636                 {
 637                 if (bitdepth == 8)
 638                     rowp[p] = outval;
 639                 else {
 640                     int rowpos = p * 2;
 641                     rowp[rowpos] = (outval >> 8) & 0xFF;
 642                     rowp[rowpos+1] = (outval & 0xFF);
 643                 }
 644             }
 645         }
 646       for (ir = 0; ir < nycell; ir++)
 647         {
 648           png_write_rows (png_ptr, &row_pointer, 1);
 649         }
 650     }
 651   sys_free(rowp, "if2img");
 652
 653   png_write_end(png_ptr, info_ptr);
 654   png_destroy_write_struct(&png_ptr, &info_ptr);
 655
 656   fclose(fp);
 657 }
 658 #endif
 659
 660 #ifdef HAVE_GD
 661 #include "gd.h"
 662 #define N_GRAYSCALE 256
 663 #endif
 664
 665 void
 666 sdf2d_to_gif (ImageFile& im, char *outfile, int nxcell, int nycell, double densmin, double densmax)
 667 {
 668 #ifdef HAVE_GD
 669   gdImagePtr gif;
 670   FILE *out;
 671   int gs_indices[N_GRAYSCALE];
 672   int i;
 673   int irow;
 674   int lastrow;
 675   unsigned char *rowp;
 676   int nx = im.nx();
 677   int ny = im.ny();
 678   ImageFileArray v = im.getArray();
 679
 680   rowp = (unsigned char *) sys_alloc(nx * nxcell,"sdf2d_to_img");
 681   if (rowp == NULL)
 682     return;
 683
 684   gif = gdImageCreate(nx * nxcell, ny * nycell);
 685   for (i = 0; i < N_GRAYSCALE; i++)
 686     gs_indices[i] = gdImageColorAllocate(gif, i, i, i);
 687
 688   lastrow = ny * nycell - 1;
 689   for (irow = 0; irow < ny; irow++)
 690     {
 691       int icol, ir;
 692       int rpos = irow * nycell;
 693       for (ir = rpos; ir < rpos + nycell; ir++)
 694         {
 695           for (icol = 0; icol < nx; icol++)
 696             {
 697               double dens;
 698               int ic;
 699               int cpos = icol * nxcell;
 700               dens = (v[icol][irow] - densmin) / (densmax - densmin);
 701               if (dens < 0)
 702                 dens = 0;
 703               else if (dens > 1)
 704                 dens = 1;
 705               for (ic = cpos; ic < cpos + nxcell; ic++)
 706                 {
 707                   rowp[ic] = (unsigned int) (dens * (double) (N_GRAYSCALE - 1));
 708                   gdImageSetPixel(gif, ic, lastrow - ir, gs_indices[rowp[ic]]);
 709                 }
 710             }
 711         }
 712     }
 713   sys_free(rowp, "if2img");
 714
 715   if ((out = fopen(outfile,"w")) == NULL)
 716     {
 717       fprintf(stderr,"Error opening output file %s for writing\n", outfile);
 718       exit(1);
 719     }
 720   gdImageGif(gif,out);
 721   fclose(out);
 722   gdImageDestroy(gif);
 723 #else
 724   fprintf(stderr, "This version does not support GIF");
 725 #endif
 726 }
 727
 728 #ifndef NO_MAIN
 729 int
 730 main (int argc, char *const argv[])
 731 {
 732   return (if2img_main(argc, argv));
 733 }
 734 #endif