home *** CD-ROM | disk | FTP | other *** search
/ rtsi.com / 2014.01.www.rtsi.com.tar / www.rtsi.com / OS9 / OSK / NETWORK / netpbm_src.lzh / NETPBM / PNM / pnmconvol.c.orig < prev    next >
Text File  |  1996-11-18  |  55KB  |  1,900 lines
  1. /* pnmconvol.c - general MxN convolution on a portable anymap
  2. **
  3. ** Version 2.0   November 26, 1994
  4. **
  5. ** Major rewriting by Mike Burns
  6. ** Copyright (C) 1994 by Mike Burns (burns@chem.psu.edu)
  7. **
  8. ** Copyright (C) 1989, 1991 by Jef Poskanzer.
  9. **
  10. ** Permission to use, copy, modify, and distribute this software and its
  11. ** documentation for any purpose and without fee is hereby granted, provided
  12. ** that the above copyright notice appear in all copies and that both that
  13. ** copyright notice and this permission notice appear in supporting
  14. ** documentation.  This software is provided "as is" without express or
  15. ** implied warranty.
  16. */
  17.  
  18. /* Version 2.0 Changes
  19. **
  20. ** Reduce run time by general optimizations and handling special cases of
  21. ** convolution matrices.  Program automatically determines if convolution 
  22. ** matrix is one of the types it can make use of so no extra command line
  23. ** arguments are necessary.
  24. **
  25. ** Examples of convolution matrices for the special cases are
  26. **
  27. **    Mean       Horizontal    Vertical
  28. **    x x x        x x x        x y z
  29. **    x x x        y y y        x y z
  30. **    x x x        z z z        x y z
  31. **
  32. ** I don't know if the horizontal and vertical ones are of much use, but
  33. ** after working on the mean convolution, it gave me ideas for the other two.
  34. **
  35. ** Some other compiler dependent optimizations
  36. ** -------------------------------------------
  37. ** Created separate functions as code was getting too large to put keep both
  38. ** PGM and PPM cases in same function and also because SWITCH statement in
  39. ** inner loop can take progressively more time the larger the size of the 
  40. ** convolution matrix.  GCC is affected this way.
  41. **
  42. ** Removed use of MOD (%) operator from innermost loop by modifying manner in
  43. ** which the current xelbuf[] is chosen.
  44. ** 
  45. */
  46.  
  47. #include "pnm.h"
  48.  
  49. /* Macros to verify that r,g,b values are within proper range */
  50.  
  51. #define CHECK_GRAY \
  52.     if ( tempgsum < 0L ) g = 0; \
  53.     else if ( tempgsum > maxval ) g = maxval; \
  54.     else g = tempgsum;
  55.  
  56. #define CHECK_RED \
  57.     if ( temprsum < 0L ) r = 0; \
  58.     else if ( temprsum > maxval ) r = maxval; \
  59.     else r = temprsum;
  60.  
  61. #define CHECK_GREEN \
  62.     if ( tempgsum < 0L ) g = 0; \
  63.     else if ( tempgsum > maxval ) g = maxval; \
  64.     else g = tempgsum;
  65.  
  66. #define CHECK_BLUE \
  67.     if ( tempbsum < 0L ) b = 0; \
  68.     else if ( tempbsum > maxval ) b = maxval; \
  69.     else b = tempbsum;
  70.  
  71.  
  72. static int check_convolve_type ARGS((xel **cxels));
  73. static void pgm_general_convolve ARGS((void));
  74. static void ppm_general_convolve ARGS((void));
  75. static void pgm_mean_convolve ARGS((void));
  76. static void ppm_mean_convolve ARGS((void));
  77. static void pgm_vertical_convolve ARGS((void));
  78. static void ppm_vertical_convolve ARGS((void));
  79. static void pgm_horizontal_convolve ARGS((void));
  80. static void ppm_horizontal_convolve ARGS((void));
  81.  
  82. #define TRUE    1
  83. #define FALSE    0
  84.  
  85. #define GENERAL_CONVOLVE    0
  86. #define MEAN_CONVOLVE        1
  87. #define HORIZONTAL_CONVOLVE    2
  88. #define VERTICAL_CONVOLVE    3
  89.  
  90. static FILE* ifp;
  91. static float** rweights;
  92. static float** gweights;
  93. static float** bweights;
  94. static int crows, ccols, ccolso2, crowso2;
  95. static int cols, rows;
  96. static xelval maxval;
  97. static int format, newformat;
  98. static float rmeanweight, gmeanweight, bmeanweight;
  99.  
  100. int
  101. main( argc, argv )
  102.     int argc;
  103.     char* argv[];
  104.     {
  105.     FILE* cifp;
  106.     xel** cxels;
  107.     int argn, cformat;
  108.     int crow, row;
  109.     register int ccol, col;
  110.     xelval cmaxval;
  111.     xelval g;
  112.     float gsum;
  113.     float rsum, bsum;
  114.     char* usage = "<convolutionfile> [pnmfile]";
  115.     int convolve_type;
  116.  
  117.     pnm_init( &argc, argv );
  118.  
  119.     argn = 1;
  120.  
  121.     if ( argn == argc )
  122.     pm_usage( usage );
  123.     cifp = pm_openr( argv[argn] );
  124.     ++argn;
  125.  
  126.     if ( argn != argc )
  127.     {
  128.     ifp = pm_openr( argv[argn] );
  129.     ++argn;
  130.     }
  131.     else
  132.     ifp = stdin;
  133.  
  134.     if ( argn != argc )
  135.     pm_usage( usage );
  136.  
  137.     pnm_pbmmaxval = PNM_MAXMAXVAL;  /* use larger value for better results */
  138.  
  139.     /* Read in the convolution matrix. */
  140.     cxels = pnm_readpnm( cifp, &ccols, &crows, &cmaxval, &cformat );
  141.     pm_close( cifp );
  142.     if ( ccols % 2 != 1 || crows % 2 != 1 )
  143.     pm_error(
  144.      "the convolution matrix must have an odd number of rows and columns" );
  145.     ccolso2 = ccols / 2;
  146.     crowso2 = crows / 2;
  147.  
  148.     pnm_readpnminit( ifp, &cols, &rows, &maxval, &format );
  149.     if ( cols < ccols || rows < crows )
  150.     pm_error(
  151.         "the image is smaller than the convolution matrix" );
  152.  
  153.     newformat = max( PNM_FORMAT_TYPE(cformat), PNM_FORMAT_TYPE(format) );
  154.     if ( PNM_FORMAT_TYPE(cformat) != newformat )
  155.     pnm_promoteformat( cxels, ccols, crows, cmaxval, cformat, cmaxval, newformat );
  156.     if ( PNM_FORMAT_TYPE(format) != newformat )
  157.         {
  158.         switch ( PNM_FORMAT_TYPE(newformat) )
  159.             {
  160.             case PPM_TYPE:
  161.             if ( PNM_FORMAT_TYPE(format) != newformat )
  162.                 pm_message( "promoting to PPM" );
  163.             break;
  164.             case PGM_TYPE:
  165.             if ( PNM_FORMAT_TYPE(format) != newformat )
  166.                 pm_message( "promoting to PGM" );
  167.             break;
  168.             }
  169.         }
  170.  
  171.  
  172.     /* Set up the normalized weights. */
  173.     rweights = (float**) pm_allocarray( ccols, crows, sizeof(float) );
  174.     gweights = (float**) pm_allocarray( ccols, crows, sizeof(float) );
  175.     bweights = (float**) pm_allocarray( ccols, crows, sizeof(float) );
  176.     rsum = gsum = bsum = 0;
  177.     for ( crow = 0; crow < crows; ++crow )
  178.     for ( ccol = 0; ccol < ccols; ++ccol )
  179.         {
  180.         switch ( PNM_FORMAT_TYPE(format) )
  181.         {
  182.         case PPM_TYPE:
  183.         rsum += rweights[crow][ccol] =
  184.             ( PPM_GETR(cxels[crow][ccol]) * 2.0 / cmaxval - 1.0 );
  185.         gsum += gweights[crow][ccol] =
  186.             ( PPM_GETG(cxels[crow][ccol]) * 2.0 / cmaxval - 1.0 );
  187.         bsum += bweights[crow][ccol] =
  188.             ( PPM_GETB(cxels[crow][ccol]) * 2.0 / cmaxval - 1.0 );
  189.         break;
  190.  
  191.         default:
  192.         gsum += gweights[crow][ccol] =
  193.             ( PNM_GET1(cxels[crow][ccol]) * 2.0 / cmaxval - 1.0 );
  194.         break;
  195.         }
  196.         }
  197.  
  198.     /* For mean_convolve routines.  All weights of a single color are the same
  199.     ** so just grab any one of them.
  200.     */
  201.     rmeanweight = rweights[0][0];
  202.     gmeanweight = gweights[0][0];
  203.     bmeanweight = bweights[0][0];
  204.  
  205.     switch ( PNM_FORMAT_TYPE(format) )
  206.     {
  207.     case PPM_TYPE:
  208.     if ( rsum < 0.9 || rsum > 1.1 || gsum < 0.9 || gsum > 1.1 ||
  209.          bsum < 0.9 || bsum > 1.1 )
  210.         pm_message(
  211.         "WARNING - this convolution matrix is biased" );
  212.     break;
  213.  
  214.     default:
  215.     if ( gsum < 0.9 || gsum > 1.1 )
  216.         pm_message(
  217.          "WARNING - this convolution matrix is biased" );
  218.     break;
  219.     }
  220.  
  221.     /* Handle certain special cases when runtime can be improved. */
  222.  
  223.     convolve_type = check_convolve_type(cxels);
  224.  
  225.     switch ( PNM_FORMAT_TYPE(format) )
  226.     {
  227.     case PPM_TYPE:
  228.     switch (convolve_type)
  229.         {
  230.         case MEAN_CONVOLVE:
  231.         ppm_mean_convolve();
  232.         break;
  233.  
  234.         case HORIZONTAL_CONVOLVE:
  235.         ppm_horizontal_convolve();
  236.         break;
  237.  
  238.         case VERTICAL_CONVOLVE:
  239.         ppm_vertical_convolve();
  240.         break;
  241.  
  242.         default: /* GENERAL_CONVOLVE */
  243.         ppm_general_convolve();
  244.         break;  
  245.         }
  246.     break;
  247.  
  248.     default:
  249.     switch (convolve_type)
  250.         {
  251.         case MEAN_CONVOLVE:
  252.         pgm_mean_convolve();
  253.         break;
  254.  
  255.         case HORIZONTAL_CONVOLVE:
  256.         pgm_horizontal_convolve();
  257.         break;
  258.  
  259.         case VERTICAL_CONVOLVE:
  260.         pgm_vertical_convolve();
  261.         break;
  262.  
  263.         default: /* GENERAL_CONVOLVE */
  264.         pgm_general_convolve();
  265.         break;  
  266.         }
  267.     break;
  268.     } /* PNM_TYPE */
  269.  
  270.     pm_close( ifp );
  271.     exit( 0 );
  272.     }
  273.  
  274.  
  275.  
  276. /* check_convolve_type
  277. **
  278. ** Determine if the convolution matrix is one of the special cases that
  279. ** can be processed faster than the general form.
  280. **
  281. ** Does not check for the case where one of the PPM colors can have differing 
  282. ** types.  Only handles cases where all PPM's are of the same special case.
  283. */
  284.  
  285. static int
  286. check_convolve_type (cxels)
  287.     xel **cxels;
  288.     {
  289.     int convolve_type;
  290.     int horizontal, vertical;
  291.     int tempcxel, rtempcxel, gtempcxel, btempcxel;
  292.     int crow, ccol;
  293.  
  294.     switch ( PNM_FORMAT_TYPE(format) )
  295.     {
  296.     case PPM_TYPE:
  297.     horizontal = TRUE;
  298.     crow = 0;
  299.     while ( horizontal && (crow < crows) )
  300.         {
  301.         ccol = 1;
  302.         rtempcxel = PPM_GETR(cxels[crow][0]);
  303.         gtempcxel = PPM_GETG(cxels[crow][0]);
  304.         btempcxel = PPM_GETB(cxels[crow][0]);
  305.         while ( horizontal && (ccol < ccols) )
  306.         {
  307.         if (( PPM_GETR(cxels[crow][ccol]) != rtempcxel ) |
  308.             ( PPM_GETG(cxels[crow][ccol]) != gtempcxel ) |
  309.             ( PPM_GETB(cxels[crow][ccol]) != btempcxel )) 
  310.             horizontal = FALSE;
  311.         ++ccol;
  312.         }
  313.         ++crow;
  314.         }
  315.  
  316.     vertical = TRUE;
  317.     ccol = 0;
  318.     while ( vertical && (ccol < ccols) )
  319.         {
  320.         crow = 1;
  321.         rtempcxel = PPM_GETR(cxels[0][ccol]);
  322.         gtempcxel = PPM_GETG(cxels[0][ccol]);
  323.         btempcxel = PPM_GETB(cxels[0][ccol]);
  324.         while ( vertical && (crow < crows) )
  325.         {
  326.         if (( PPM_GETR(cxels[crow][ccol]) != rtempcxel ) |
  327.             ( PPM_GETG(cxels[crow][ccol]) != gtempcxel ) |
  328.             ( PPM_GETB(cxels[crow][ccol]) != btempcxel ))
  329.             vertical = FALSE;
  330.         ++crow;
  331.         }
  332.         ++ccol;
  333.         }
  334.     break;
  335.  
  336.     default:
  337.     horizontal = TRUE;
  338.     crow = 0;
  339.     while ( horizontal && (crow < crows) )
  340.         {
  341.         ccol = 1;
  342.         tempcxel = PNM_GET1(cxels[crow][0]);
  343.         while ( horizontal && (ccol < ccols) )
  344.         {
  345.         if ( PNM_GET1(cxels[crow][ccol]) != tempcxel )
  346.             horizontal = FALSE;
  347.         ++ccol;
  348.         }
  349.         ++crow;
  350.         }
  351.  
  352.     vertical = TRUE;
  353.     ccol = 0;
  354.     while ( vertical && (ccol < ccols) )
  355.         {
  356.         crow = 1;
  357.         tempcxel = PNM_GET1(cxels[0][ccol]);
  358.         while ( vertical && (crow < crows) )
  359.         {
  360.         if ( PNM_GET1(cxels[crow][ccol]) != tempcxel )
  361.             vertical = FALSE;
  362.         ++crow;
  363.         }
  364.         ++ccol;
  365.         }
  366.     break;
  367.     }
  368.  
  369.     /* Which type do we have? */
  370.     if ( horizontal && vertical )
  371.     convolve_type = MEAN_CONVOLVE;
  372.     else if ( horizontal )
  373.     convolve_type = HORIZONTAL_CONVOLVE;
  374.     else if ( vertical )
  375.     convolve_type = VERTICAL_CONVOLVE;
  376.     else
  377.     convolve_type = GENERAL_CONVOLVE;
  378.  
  379.     return (convolve_type);
  380.     }
  381.  
  382.  
  383.  
  384.  
  385. /* General PGM Convolution
  386. **
  387. ** No useful redundancy in convolution matrix.
  388. */
  389.  
  390. static void
  391. pgm_general_convolve()
  392.     {
  393.     register int ccol, col;
  394.     xel** xelbuf;
  395.     xel* outputrow;
  396.     xel x, y;
  397.     xelval g;
  398.     int row, crow;
  399.     float gsum;
  400.     xel **rowptr, *temprptr;
  401.     int toprow, temprow;
  402.     int i, irow;
  403.     int leftcol;
  404.     long tempgsum;
  405.  
  406.     /* Allocate space for one convolution-matrix's worth of rows, plus
  407.     ** a row output buffer. */
  408.     xelbuf = pnm_allocarray( cols, crows );
  409.     outputrow = pnm_allocrow( cols );
  410.  
  411.     /* Allocate array of pointers to xelbuf */
  412.     rowptr = (xel **) pnm_allocarray( 1, crows );
  413.  
  414.     pnm_writepnminit( stdout, cols, rows, maxval, newformat, 0 );
  415.  
  416.     /* Read in one convolution-matrix's worth of image, less one row. */
  417.     for ( row = 0; row < crows - 1; ++row )
  418.     {
  419.     pnm_readpnmrow( ifp, xelbuf[row], cols, maxval, format );
  420.     if ( PNM_FORMAT_TYPE(format) != newformat )
  421.         pnm_promoteformatrow(
  422.         xelbuf[row], cols, maxval, format, maxval, newformat );
  423.     /* Write out just the part we're not going to convolve. */
  424.     if ( row < crowso2 )
  425.         pnm_writepnmrow( stdout, xelbuf[row], cols, maxval, newformat, 0 );
  426.     }
  427.  
  428.     /* Now the rest of the image - read in the row at the end of
  429.     ** xelbuf, and convolve and write out the row in the middle.
  430.     */
  431.     for ( ; row < rows; ++row )
  432.     {
  433.     toprow = row + 1;
  434.     temprow = row % crows;
  435.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  436.     if ( PNM_FORMAT_TYPE(format) != newformat )
  437.         pnm_promoteformatrow(
  438.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  439.  
  440.     /* Arrange rowptr to eliminate the use of mod function to determine
  441.     ** which row of xelbuf is 0...crows.  Mod function can be very costly.
  442.     */
  443.     temprow = toprow % crows;
  444.     i = 0;
  445.     for (irow = temprow; irow < crows; ++i, ++irow)
  446.         rowptr[i] = xelbuf[irow];
  447.     for (irow = 0; irow < temprow; ++irow, ++i)
  448.         rowptr[i] = xelbuf[irow];
  449.  
  450.     for ( col = 0; col < cols; ++col )
  451.         {
  452.         if ( col < ccolso2 || col >= cols - ccolso2 )
  453.         outputrow[col] = rowptr[crowso2][col];
  454.         else
  455.         {
  456.         leftcol = col - ccolso2;
  457.         gsum = 0.0;
  458.         for ( crow = 0; crow < crows; ++crow )
  459.             {
  460.             temprptr = rowptr[crow] + leftcol;
  461.             for ( ccol = 0; ccol < ccols; ++ccol )
  462.             gsum += PNM_GET1( *(temprptr + ccol) )
  463.                 * gweights[crow][ccol];
  464.             }
  465.             tempgsum = gsum + 0.5;
  466.             CHECK_GRAY;
  467.             PNM_ASSIGN1( outputrow[col], g );
  468.         }
  469.         }
  470.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  471.     }
  472.  
  473.     /* Now write out the remaining unconvolved rows in xelbuf. */
  474.     for ( irow = crowso2 + 1; irow < crows; ++irow )
  475.     pnm_writepnmrow(
  476.             stdout, rowptr[irow], cols, maxval, newformat, 0 );
  477.  
  478.     pm_close( stdout );
  479.     }
  480.  
  481.  
  482. /* PGM Mean Convolution
  483. **
  484. ** This is the common case where you just want the target pixel replaced with
  485. ** the average value of its neighbors.  This can work much faster than the
  486. ** general case because you can reduce the number of floating point operations
  487. ** that are required since all the weights are the same.  You will only need
  488. ** to multiply by the weight once, not for every pixel in the convolution
  489. ** matrix.
  490. **
  491. ** This algorithm works by creating sums for each column of crows height for
  492. ** the whole width of the image.  Then add ccols column sums together to obtain
  493. ** the total sum of the neighbors and multiply that sum by the weight.  As you
  494. ** move right to left to calculate the next pixel, take the total sum you just
  495. ** generated, add in the value of the next column and subtract the value of the
  496. ** leftmost column.  Multiply that by the weight and that's it.  As you move
  497. ** down a row, calculate new column sums by using previous sum for that column
  498. ** and adding in pixel on current row and subtracting pixel in top row.
  499. **
  500. */
  501.  
  502.  
  503. static void
  504. pgm_mean_convolve()
  505.     {
  506.     register int ccol, col;
  507.     xel** xelbuf;
  508.     xel* outputrow;
  509.     xelval g;
  510.     int row, crow;
  511.     xel **rowptr, *temprptr;
  512.     int leftcol;
  513.     int i, irow;
  514.     int toprow, temprow;
  515.     int subrow, addrow;
  516.     int subcol, addcol;
  517.     long gisum;
  518.     int tempcol, crowsp1;
  519.     long tempgsum;
  520.     long *gcolumnsum;
  521.  
  522.     /* Allocate space for one convolution-matrix's worth of rows, plus
  523.     ** a row output buffer.  MEAN uses an extra row. */
  524.     xelbuf = pnm_allocarray( cols, crows + 1 );
  525.     outputrow = pnm_allocrow( cols );
  526.  
  527.     /* Allocate array of pointers to xelbuf. MEAN uses an extra row. */
  528.     rowptr = (xel **) pnm_allocarray( 1, crows + 1);
  529.  
  530.     /* Allocate space for intermediate column sums */
  531.     gcolumnsum = (long *) pm_allocrow( cols, sizeof(long) );
  532.     for ( col = 0; col < cols; ++col )
  533.     gcolumnsum[col] = 0L;
  534.  
  535.     pnm_writepnminit( stdout, cols, rows, maxval, newformat, 0 );
  536.  
  537.     /* Read in one convolution-matrix's worth of image, less one row. */
  538.     for ( row = 0; row < crows - 1; ++row )
  539.     {
  540.     pnm_readpnmrow( ifp, xelbuf[row], cols, maxval, format );
  541.     if ( PNM_FORMAT_TYPE(format) != newformat )
  542.         pnm_promoteformatrow(
  543.         xelbuf[row], cols, maxval, format, maxval, newformat );
  544.     /* Write out just the part we're not going to convolve. */
  545.     if ( row < crowso2 )
  546.         pnm_writepnmrow( stdout, xelbuf[row], cols, maxval, newformat, 0 );
  547.     }
  548.  
  549.     /* Do first real row only */
  550.     subrow = crows;
  551.     addrow = crows - 1;
  552.     toprow = row + 1;
  553.     temprow = row % crows;
  554.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  555.     if ( PNM_FORMAT_TYPE(format) != newformat )
  556.     pnm_promoteformatrow(
  557.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  558.  
  559.     temprow = toprow % crows;
  560.     i = 0;
  561.     for (irow = temprow; irow < crows; ++i, ++irow)
  562.     rowptr[i] = xelbuf[irow];
  563.     for (irow = 0; irow < temprow; ++irow, ++i)
  564.     rowptr[i] = xelbuf[irow];
  565.  
  566.     gisum = 0L;
  567.     for ( col = 0; col < cols; ++col )
  568.     {
  569.     if ( col < ccolso2 || col >= cols - ccolso2 )
  570.         outputrow[col] = rowptr[crowso2][col];
  571.     else if ( col == ccolso2 )
  572.         {
  573.         leftcol = col - ccolso2;
  574.         for ( crow = 0; crow < crows; ++crow )
  575.         {
  576.         temprptr = rowptr[crow] + leftcol;
  577.         for ( ccol = 0; ccol < ccols; ++ccol )
  578.             gcolumnsum[leftcol + ccol] += 
  579.             PNM_GET1( *(temprptr + ccol) );
  580.         }
  581.         for ( ccol = 0; ccol < ccols; ++ccol)
  582.         gisum += gcolumnsum[leftcol + ccol];
  583.         tempgsum = (float) gisum * gmeanweight + 0.5;
  584.         CHECK_GRAY;
  585.         PNM_ASSIGN1( outputrow[col], g );
  586.         }
  587.     else
  588.         {
  589.         /* Column numbers to subtract or add to isum */
  590.         subcol = col - ccolso2 - 1;
  591.         addcol = col + ccolso2;  
  592.         for ( crow = 0; crow < crows; ++crow )
  593.         gcolumnsum[addcol] += PNM_GET1( rowptr[crow][addcol] );
  594.         gisum = gisum - gcolumnsum[subcol] + gcolumnsum[addcol];
  595.         tempgsum = (float) gisum * gmeanweight + 0.5;
  596.         CHECK_GRAY;
  597.         PNM_ASSIGN1( outputrow[col], g );
  598.         }
  599.     }
  600.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  601.  
  602.     ++row;
  603.     /* For all subsequent rows do it this way as the columnsums have been
  604.     ** generated.  Now we can use them to reduce further calculations.
  605.     */
  606.     crowsp1 = crows + 1;
  607.     for ( ; row < rows; ++row )
  608.     {
  609.     toprow = row + 1;
  610.     temprow = row % (crows + 1);
  611.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  612.     if ( PNM_FORMAT_TYPE(format) != newformat )
  613.         pnm_promoteformatrow(
  614.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  615.  
  616.     /* This rearrangement using crows+1 rowptrs and xelbufs will cause
  617.     ** rowptr[0..crows-1] to always hold active xelbufs and for 
  618.     ** rowptr[crows] to always hold the oldest (top most) xelbuf.
  619.     */
  620.     temprow = (toprow + 1) % crowsp1;
  621.     i = 0;
  622.     for (irow = temprow; irow < crowsp1; ++i, ++irow)
  623.         rowptr[i] = xelbuf[irow];
  624.     for (irow = 0; irow < temprow; ++irow, ++i)
  625.         rowptr[i] = xelbuf[irow];
  626.  
  627.     gisum = 0L;
  628.     for ( col = 0; col < cols; ++col )
  629.         {
  630.         if ( col < ccolso2 || col >= cols - ccolso2 )
  631.         outputrow[col] = rowptr[crowso2][col];
  632.         else if ( col == ccolso2 )
  633.         {
  634.         leftcol = col - ccolso2;
  635.         for ( ccol = 0; ccol < ccols; ++ccol )
  636.             {
  637.             tempcol = leftcol + ccol;
  638.             gcolumnsum[tempcol] = gcolumnsum[tempcol]
  639.             - PNM_GET1( rowptr[subrow][ccol] )
  640.             + PNM_GET1( rowptr[addrow][ccol] );
  641.             gisum += gcolumnsum[tempcol];
  642.             }
  643.         tempgsum = (float) gisum * gmeanweight + 0.5;
  644.         CHECK_GRAY;
  645.         PNM_ASSIGN1( outputrow[col], g );
  646.         }
  647.         else
  648.         {
  649.         /* Column numbers to subtract or add to isum */
  650.         subcol = col - ccolso2 - 1;
  651.         addcol = col + ccolso2;  
  652.         gcolumnsum[addcol] = gcolumnsum[addcol]
  653.             - PNM_GET1( rowptr[subrow][addcol] )
  654.             + PNM_GET1( rowptr[addrow][addcol] );
  655.         gisum = gisum - gcolumnsum[subcol] + gcolumnsum[addcol];
  656.         tempgsum = (float) gisum * gmeanweight + 0.5;
  657.         CHECK_GRAY;
  658.         PNM_ASSIGN1( outputrow[col], g );
  659.         }
  660.         }
  661.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  662.     }
  663.  
  664.     /* Now write out the remaining unconvolved rows in xelbuf. */
  665.     for ( irow = crowso2 + 1; irow < crows; ++irow )
  666.     pnm_writepnmrow(
  667.             stdout, rowptr[irow], cols, maxval, newformat, 0 );
  668.  
  669.     pm_close( stdout );
  670.     }
  671.  
  672.  
  673. /* PGM Horizontal Convolution
  674. **
  675. ** Similar idea to using columnsums of the Mean and Vertical convolution,
  676. ** but uses temporary sums of row values.  Need to multiply by weights crows
  677. ** number of times.  Each time a new line is started, must recalculate the
  678. ** initials rowsums for the newest row only.  Uses queue to still access
  679. ** previous row sums.
  680. **
  681. */
  682.  
  683. static void
  684. pgm_horizontal_convolve()
  685.     {
  686.     register int ccol, col;
  687.     xel** xelbuf;
  688.     xel* outputrow;
  689.     xel x;
  690.     xelval g;
  691.     int row, crow;
  692.     xel **rowptr, *temprptr;
  693.     int leftcol;
  694.     int i, irow;
  695.     int temprow;
  696.     int subcol, addcol;
  697.     float gsum;
  698.     int addrow, subrow;
  699.     long **growsum, **growsumptr;
  700.     int crowsp1;
  701.     long tempgsum;
  702.  
  703.     /* Allocate space for one convolution-matrix's worth of rows, plus
  704.     ** a row output buffer. */
  705.     xelbuf = pnm_allocarray( cols, crows + 1 );
  706.     outputrow = pnm_allocrow( cols );
  707.  
  708.     /* Allocate array of pointers to xelbuf */
  709.     rowptr = (xel **) pnm_allocarray( 1, crows + 1);
  710.  
  711.     /* Allocate intermediate row sums.  HORIZONTAL uses an extra row. */
  712.     /* crows current rows and 1 extra for newest added row.           */
  713.     growsum = (long **) pm_allocarray( cols, crows + 1, sizeof(float) );
  714.     growsumptr = (long **) pnm_allocarray( 1, crows + 1);
  715.  
  716.     pnm_writepnminit( stdout, cols, rows, maxval, newformat, 0 );
  717.  
  718.     /* Read in one convolution-matrix's worth of image, less one row. */
  719.     for ( row = 0; row < crows - 1; ++row )
  720.     {
  721.     pnm_readpnmrow( ifp, xelbuf[row], cols, maxval, format );
  722.     if ( PNM_FORMAT_TYPE(format) != newformat )
  723.         pnm_promoteformatrow(
  724.         xelbuf[row], cols, maxval, format, maxval, newformat );
  725.     /* Write out just the part we're not going to convolve. */
  726.     if ( row < crowso2 )
  727.         pnm_writepnmrow( stdout, xelbuf[row], cols, maxval, newformat, 0 );
  728.     }
  729.  
  730.     /* First row only */
  731.     temprow = row % crows;
  732.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  733.     if ( PNM_FORMAT_TYPE(format) != newformat )
  734.     pnm_promoteformatrow(
  735.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  736.  
  737.     temprow = (row + 1) % crows;
  738.     i = 0;
  739.     for (irow = temprow; irow < crows; ++i, ++irow)
  740.     rowptr[i] = xelbuf[irow];
  741.     for (irow = 0; irow < temprow; ++irow, ++i)
  742.     rowptr[i] = xelbuf[irow];
  743.  
  744.     for ( crow = 0; crow < crows; ++crow )
  745.     growsumptr[crow] = growsum[crow];
  746.  
  747.     for ( col = 0; col < cols; ++col )
  748.     {
  749.     if ( col < ccolso2 || col >= cols - ccolso2 )
  750.         outputrow[col] = rowptr[crowso2][col];
  751.     else if ( col == ccolso2 )
  752.         {
  753.         leftcol = col - ccolso2;
  754.         gsum = 0.0;
  755.         for ( crow = 0; crow < crows; ++crow )
  756.         {
  757.         temprptr = rowptr[crow] + leftcol;
  758.         growsumptr[crow][leftcol] = 0L;
  759.         for ( ccol = 0; ccol < ccols; ++ccol )
  760.             growsumptr[crow][leftcol] += 
  761.                 PNM_GET1( *(temprptr + ccol) );
  762.         gsum += growsumptr[crow][leftcol] * gweights[crow][0];
  763.         }
  764.         tempgsum = gsum + 0.5;
  765.         CHECK_GRAY;
  766.         PNM_ASSIGN1( outputrow[col], g );
  767.         }
  768.     else
  769.         {
  770.         gsum = 0.0;
  771.         leftcol = col - ccolso2;
  772.         subcol = col - ccolso2 - 1;
  773.         addcol = col + ccolso2;
  774.         for ( crow = 0; crow < crows; ++crow )
  775.         {
  776.         growsumptr[crow][leftcol] = growsumptr[crow][subcol]
  777.             - PNM_GET1( rowptr[crow][subcol] )
  778.             + PNM_GET1( rowptr[crow][addcol] );
  779.         gsum += growsumptr[crow][leftcol] * gweights[crow][0];
  780.         }
  781.         tempgsum = gsum + 0.5;
  782.         CHECK_GRAY;
  783.         PNM_ASSIGN1( outputrow[col], g );
  784.         }
  785.         }
  786.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  787.  
  788.  
  789.     /* For all subsequent rows */
  790.  
  791.     subrow = crows;
  792.     addrow = crows - 1;
  793.     crowsp1 = crows + 1;
  794.     ++row;
  795.     for ( ; row < rows; ++row )
  796.     {
  797.     temprow = row % crowsp1;
  798.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  799.     if ( PNM_FORMAT_TYPE(format) != newformat )
  800.         pnm_promoteformatrow(
  801.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  802.  
  803.     temprow = (row + 2) % crowsp1;
  804.     i = 0;
  805.     for (irow = temprow; irow < crowsp1; ++i, ++irow)
  806.         {
  807.         rowptr[i] = xelbuf[irow];
  808.         growsumptr[i] = growsum[irow];
  809.         }
  810.     for (irow = 0; irow < temprow; ++irow, ++i)
  811.         {
  812.         rowptr[i] = xelbuf[irow];
  813.         growsumptr[i] = growsum[irow];
  814.         }
  815.  
  816.     for ( col = 0; col < cols; ++col )
  817.         {
  818.         if ( col < ccolso2 || col >= cols - ccolso2 )
  819.         outputrow[col] = rowptr[crowso2][col];
  820.         else if ( col == ccolso2 )
  821.         {
  822.         gsum = 0.0;
  823.         leftcol = col - ccolso2;
  824.         growsumptr[addrow][leftcol] = 0L;
  825.         for ( ccol = 0; ccol < ccols; ++ccol )
  826.             growsumptr[addrow][leftcol] += 
  827.             PNM_GET1( rowptr[addrow][leftcol + ccol] );
  828.         for ( crow = 0; crow < crows; ++crow )
  829.             gsum += growsumptr[crow][leftcol] * gweights[crow][0];
  830.         tempgsum = gsum + 0.5;
  831.         CHECK_GRAY;
  832.         PNM_ASSIGN1( outputrow[col], g );
  833.         }
  834.         else
  835.         {
  836.         gsum = 0.0;
  837.         leftcol = col - ccolso2;
  838.         subcol = col - ccolso2 - 1;
  839.         addcol = col + ccolso2;  
  840.         growsumptr[addrow][leftcol] = growsumptr[addrow][subcol]
  841.             - PNM_GET1( rowptr[addrow][subcol] )
  842.             + PNM_GET1( rowptr[addrow][addcol] );
  843.         for ( crow = 0; crow < crows; ++crow )
  844.             gsum += growsumptr[crow][leftcol] * gweights[crow][0];
  845.         tempgsum = gsum + 0.5;
  846.         CHECK_GRAY;
  847.         PNM_ASSIGN1( outputrow[col], g );
  848.         }
  849.         }
  850.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  851.     }
  852.  
  853.     /* Now write out the remaining unconvolved rows in xelbuf. */
  854.     for ( irow = crowso2 + 1; irow < crows; ++irow )
  855.     pnm_writepnmrow(
  856.             stdout, rowptr[irow], cols, maxval, newformat, 0 );
  857.  
  858.     pm_close( stdout );
  859.     }
  860.  
  861.  
  862. /* PGM Vertical Convolution
  863. **
  864. ** Uses column sums as in Mean Convolution.
  865. **
  866. */
  867.  
  868.  
  869. static void
  870. pgm_vertical_convolve()
  871.     {
  872.     register int ccol, col;
  873.     xel** xelbuf;
  874.     xel* outputrow;
  875.     xelval g;
  876.     int row, crow;
  877.     xel **rowptr, *temprptr;
  878.     int leftcol;
  879.     int i, irow;
  880.     int toprow, temprow;
  881.     int subrow, addrow;
  882.     int tempcol;
  883.     float gsum;
  884.     long *gcolumnsum;
  885.     int crowsp1;
  886.     int addcol;
  887.     long tempgsum;
  888.  
  889.     /* Allocate space for one convolution-matrix's worth of rows, plus
  890.     ** a row output buffer. VERTICAL uses an extra row. */
  891.     xelbuf = pnm_allocarray( cols, crows + 1 );
  892.     outputrow = pnm_allocrow( cols );
  893.  
  894.     /* Allocate array of pointers to xelbuf */
  895.     rowptr = (xel **) pnm_allocarray( 1, crows + 1 );
  896.  
  897.     /* Allocate space for intermediate column sums */
  898.     gcolumnsum = (long *) pm_allocrow( cols, sizeof(long) );
  899.     for ( col = 0; col < cols; ++col )
  900.     gcolumnsum[col] = 0L;
  901.  
  902.     pnm_writepnminit( stdout, cols, rows, maxval, newformat, 0 );
  903.  
  904.     /* Read in one convolution-matrix's worth of image, less one row. */
  905.     for ( row = 0; row < crows - 1; ++row )
  906.     {
  907.     pnm_readpnmrow( ifp, xelbuf[row], cols, maxval, format );
  908.     if ( PNM_FORMAT_TYPE(format) != newformat )
  909.         pnm_promoteformatrow(
  910.         xelbuf[row], cols, maxval, format, maxval, newformat );
  911.     /* Write out just the part we're not going to convolve. */
  912.     if ( row < crowso2 )
  913.         pnm_writepnmrow( stdout, xelbuf[row], cols, maxval, newformat, 0 );
  914.     }
  915.  
  916.     /* Now the rest of the image - read in the row at the end of
  917.     ** xelbuf, and convolve and write out the row in the middle.
  918.     */
  919.     /* For first row only */
  920.  
  921.     toprow = row + 1;
  922.     temprow = row % crows;
  923.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  924.     if ( PNM_FORMAT_TYPE(format) != newformat )
  925.     pnm_promoteformatrow(
  926.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  927.  
  928.     /* Arrange rowptr to eliminate the use of mod function to determine
  929.     ** which row of xelbuf is 0...crows.  Mod function can be very costly.
  930.     */
  931.     temprow = toprow % crows;
  932.     i = 0;
  933.     for (irow = temprow; irow < crows; ++i, ++irow)
  934.     rowptr[i] = xelbuf[irow];
  935.     for (irow = 0; irow < temprow; ++irow, ++i)
  936.     rowptr[i] = xelbuf[irow];
  937.  
  938.     for ( col = 0; col < cols; ++col )
  939.     {
  940.     if ( col < ccolso2 || col >= cols - ccolso2 )
  941.         outputrow[col] = rowptr[crowso2][col];
  942.     else if ( col == ccolso2 )
  943.         {
  944.         gsum = 0.0;
  945.         leftcol = col - ccolso2;
  946.         for ( crow = 0; crow < crows; ++crow )
  947.         {
  948.         temprptr = rowptr[crow] + leftcol;
  949.         for ( ccol = 0; ccol < ccols; ++ccol )
  950.             gcolumnsum[leftcol + ccol] += 
  951.             PNM_GET1( *(temprptr + ccol) );
  952.         }
  953.         for ( ccol = 0; ccol < ccols; ++ccol)
  954.         gsum += gcolumnsum[leftcol + ccol] * gweights[0][ccol];
  955.         tempgsum = gsum + 0.5;
  956.         CHECK_GRAY;
  957.         PNM_ASSIGN1( outputrow[col], g );
  958.         }
  959.     else
  960.         {
  961.         gsum = 0.0;
  962.         leftcol = col - ccolso2;
  963.         addcol = col + ccolso2;  
  964.         for ( crow = 0; crow < crows; ++crow )
  965.         gcolumnsum[addcol] += PNM_GET1( rowptr[crow][addcol] );
  966.         for ( ccol = 0; ccol < ccols; ++ccol )
  967.         gsum += gcolumnsum[leftcol + ccol] * gweights[0][ccol];
  968.         tempgsum = gsum + 0.5;
  969.         CHECK_GRAY;
  970.         PNM_ASSIGN1( outputrow[col], g );
  971.         }
  972.     }
  973.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  974.  
  975.     /* For all subsequent rows */
  976.     subrow = crows;
  977.     addrow = crows - 1;
  978.     crowsp1 = crows + 1;
  979.     ++row;
  980.     for ( ; row < rows; ++row )
  981.     {
  982.     toprow = row + 1;
  983.     temprow = row % (crows +1);
  984.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  985.     if ( PNM_FORMAT_TYPE(format) != newformat )
  986.         pnm_promoteformatrow(
  987.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  988.  
  989.     /* Arrange rowptr to eliminate the use of mod function to determine
  990.     ** which row of xelbuf is 0...crows.  Mod function can be very costly.
  991.     */
  992.     temprow = (toprow + 1) % crowsp1;
  993.     i = 0;
  994.     for (irow = temprow; irow < crowsp1; ++i, ++irow)
  995.         rowptr[i] = xelbuf[irow];
  996.     for (irow = 0; irow < temprow; ++irow, ++i)
  997.         rowptr[i] = xelbuf[irow];
  998.  
  999.     for ( col = 0; col < cols; ++col )
  1000.         {
  1001.         if ( col < ccolso2 || col >= cols - ccolso2 )
  1002.         outputrow[col] = rowptr[crowso2][col];
  1003.         else if ( col == ccolso2 )
  1004.         {
  1005.         gsum = 0.0;
  1006.         leftcol = col - ccolso2;
  1007.         for ( ccol = 0; ccol < ccols; ++ccol )
  1008.             {
  1009.             tempcol = leftcol + ccol;
  1010.             gcolumnsum[tempcol] = gcolumnsum[tempcol] 
  1011.             - PNM_GET1( rowptr[subrow][ccol] )
  1012.             + PNM_GET1( rowptr[addrow][ccol] );
  1013.             gsum = gsum + gcolumnsum[tempcol] * gweights[0][ccol];
  1014.             }
  1015.         tempgsum = gsum + 0.5;
  1016.         CHECK_GRAY;
  1017.         PNM_ASSIGN1( outputrow[col], g );
  1018.         }
  1019.         else
  1020.         {
  1021.         gsum = 0.0;
  1022.         leftcol = col - ccolso2;
  1023.         addcol = col + ccolso2;
  1024.         gcolumnsum[addcol] = gcolumnsum[addcol]
  1025.             - PNM_GET1( rowptr[subrow][addcol] )
  1026.             + PNM_GET1( rowptr[addrow][addcol] );
  1027.         for ( ccol = 0; ccol < ccols; ++ccol )
  1028.             gsum += gcolumnsum[leftcol + ccol] * gweights[0][ccol];
  1029.         tempgsum = gsum + 0.5;
  1030.         CHECK_GRAY;
  1031.         PNM_ASSIGN1( outputrow[col], g );
  1032.         }
  1033.         }
  1034.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  1035.     }
  1036.  
  1037.     /* Now write out the remaining unconvolved rows in xelbuf. */
  1038.     for ( irow = crowso2 + 1; irow < crows; ++irow )
  1039.     pnm_writepnmrow(
  1040.             stdout, rowptr[irow], cols, maxval, newformat, 0 );
  1041.  
  1042.     pm_close( stdout );
  1043.     }
  1044.  
  1045.  
  1046.  
  1047.  
  1048. /* PPM General Convolution Algorithm
  1049. **
  1050. ** No redundancy in convolution matrix.  Just use brute force.
  1051. ** See pgm_general_convolve() for more details.
  1052. */
  1053.  
  1054. static void
  1055. ppm_general_convolve()
  1056.     {
  1057.     register int ccol, col;
  1058.     xel** xelbuf;
  1059.     xel* outputrow;
  1060.     xel x, y;
  1061.     xelval r, g, b;
  1062.     int row, crow;
  1063.     float rsum, gsum, bsum;
  1064.     xel **rowptr, *temprptr;
  1065.     int toprow, temprow;
  1066.     int i, irow;
  1067.     int leftcol;
  1068.     long temprsum, tempgsum, tempbsum;
  1069.  
  1070.     /* Allocate space for one convolution-matrix's worth of rows, plus
  1071.     ** a row output buffer. */
  1072.     xelbuf = pnm_allocarray( cols, crows );
  1073.     outputrow = pnm_allocrow( cols );
  1074.  
  1075.     /* Allocate array of pointers to xelbuf */
  1076.     rowptr = (xel **) pnm_allocarray( 1, crows );
  1077.  
  1078.     pnm_writepnminit( stdout, cols, rows, maxval, newformat, 0 );
  1079.  
  1080.     /* Read in one convolution-matrix's worth of image, less one row. */
  1081.     for ( row = 0; row < crows - 1; ++row )
  1082.     {
  1083.     pnm_readpnmrow( ifp, xelbuf[row], cols, maxval, format );
  1084.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1085.         pnm_promoteformatrow(
  1086.         xelbuf[row], cols, maxval, format, maxval, newformat );
  1087.     /* Write out just the part we're not going to convolve. */
  1088.     if ( row < crowso2 )
  1089.         pnm_writepnmrow( stdout, xelbuf[row], cols, maxval, newformat, 0 );
  1090.     }
  1091.  
  1092.     /* Now the rest of the image - read in the row at the end of
  1093.     ** xelbuf, and convolve and write out the row in the middle.
  1094.     */
  1095.     for ( ; row < rows; ++row )
  1096.     {
  1097.     toprow = row + 1;
  1098.     temprow = row % crows;
  1099.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  1100.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1101.         pnm_promoteformatrow(
  1102.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  1103.  
  1104.     /* Arrange rowptr to eliminate the use of mod function to determine
  1105.     ** which row of xelbuf is 0...crows.  Mod function can be very costly.
  1106.     */
  1107.     temprow = toprow % crows;
  1108.     i = 0;
  1109.     for (irow = temprow; irow < crows; ++i, ++irow)
  1110.         rowptr[i] = xelbuf[irow];
  1111.     for (irow = 0; irow < temprow; ++irow, ++i)
  1112.         rowptr[i] = xelbuf[irow];
  1113.  
  1114.     for ( col = 0; col < cols; ++col )
  1115.         {
  1116.         if ( col < ccolso2 || col >= cols - ccolso2 )
  1117.         outputrow[col] = rowptr[crowso2][col];
  1118.         else
  1119.         {
  1120.         leftcol = col - ccolso2;
  1121.         rsum = gsum = bsum = 0.0;
  1122.         for ( crow = 0; crow < crows; ++crow )
  1123.             {
  1124.             temprptr = rowptr[crow] + leftcol;
  1125.             for ( ccol = 0; ccol < ccols; ++ccol )
  1126.             {
  1127.             rsum += PPM_GETR( *(temprptr + ccol) )
  1128.                 * rweights[crow][ccol];
  1129.             gsum += PPM_GETG( *(temprptr + ccol) )
  1130.                 * gweights[crow][ccol];
  1131.             bsum += PPM_GETB( *(temprptr + ccol) )
  1132.                 * bweights[crow][ccol];
  1133.             }
  1134.             }
  1135.             temprsum = rsum + 0.5;
  1136.             tempgsum = gsum + 0.5;
  1137.             tempbsum = bsum + 0.5;
  1138.             CHECK_RED;
  1139.             CHECK_GREEN;
  1140.             CHECK_BLUE;
  1141.             PPM_ASSIGN( outputrow[col], r, g, b );
  1142.         }
  1143.         }
  1144.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  1145.     }
  1146.  
  1147.     /* Now write out the remaining unconvolved rows in xelbuf. */
  1148.     for ( irow = crowso2 + 1; irow < crows; ++irow )
  1149.     pnm_writepnmrow(
  1150.             stdout, rowptr[irow], cols, maxval, newformat, 0 );
  1151.  
  1152.     pm_close( stdout );
  1153.     }
  1154.  
  1155.  
  1156. /* PPM Mean Convolution
  1157. **
  1158. ** Same as pgm_mean_convolve() but for PPM.
  1159. **
  1160. */
  1161.  
  1162. static void
  1163. ppm_mean_convolve()
  1164.     {
  1165.     register int ccol, col;
  1166.     xel** xelbuf;
  1167.     xel* outputrow;
  1168.     xelval r, g, b;
  1169.     int row, crow;
  1170.     xel **rowptr, *temprptr;
  1171.     int leftcol;
  1172.     int i, irow;
  1173.     int toprow, temprow;
  1174.     int subrow, addrow;
  1175.     int subcol, addcol;
  1176.     long risum, gisum, bisum;
  1177.     float rsum, gsum, bsum;
  1178.     long temprsum, tempgsum, tempbsum;
  1179.     int tempcol, crowsp1;
  1180.     long *rcolumnsum, *gcolumnsum, *bcolumnsum;
  1181.  
  1182.     /* Allocate space for one convolution-matrix's worth of rows, plus
  1183.     ** a row output buffer.  MEAN uses an extra row. */
  1184.     xelbuf = pnm_allocarray( cols, crows + 1 );
  1185.     outputrow = pnm_allocrow( cols );
  1186.  
  1187.     /* Allocate array of pointers to xelbuf. MEAN uses an extra row. */
  1188.     rowptr = (xel **) pnm_allocarray( 1, crows + 1);
  1189.  
  1190.     /* Allocate space for intermediate column sums */
  1191.     rcolumnsum = (long *) pm_allocrow( cols, sizeof(long) );
  1192.     gcolumnsum = (long *) pm_allocrow( cols, sizeof(long) );
  1193.     bcolumnsum = (long *) pm_allocrow( cols, sizeof(long) );
  1194.     for ( col = 0; col < cols; ++col )
  1195.     {
  1196.     rcolumnsum[col] = 0L;
  1197.     gcolumnsum[col] = 0L;
  1198.     bcolumnsum[col] = 0L;
  1199.     }
  1200.  
  1201.     pnm_writepnminit( stdout, cols, rows, maxval, newformat, 0 );
  1202.  
  1203.     /* Read in one convolution-matrix's worth of image, less one row. */
  1204.     for ( row = 0; row < crows - 1; ++row )
  1205.     {
  1206.     pnm_readpnmrow( ifp, xelbuf[row], cols, maxval, format );
  1207.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1208.         pnm_promoteformatrow(
  1209.         xelbuf[row], cols, maxval, format, maxval, newformat );
  1210.     /* Write out just the part we're not going to convolve. */
  1211.     if ( row < crowso2 )
  1212.         pnm_writepnmrow( stdout, xelbuf[row], cols, maxval, newformat, 0 );
  1213.     }
  1214.  
  1215.     /* Do first real row only */
  1216.     subrow = crows;
  1217.     addrow = crows - 1;
  1218.     toprow = row + 1;
  1219.     temprow = row % crows;
  1220.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  1221.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1222.     pnm_promoteformatrow(
  1223.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  1224.  
  1225.     temprow = toprow % crows;
  1226.     i = 0;
  1227.     for (irow = temprow; irow < crows; ++i, ++irow)
  1228.     rowptr[i] = xelbuf[irow];
  1229.     for (irow = 0; irow < temprow; ++irow, ++i)
  1230.     rowptr[i] = xelbuf[irow];
  1231.  
  1232.     risum = 0L;
  1233.     gisum = 0L;
  1234.     bisum = 0L;
  1235.     for ( col = 0; col < cols; ++col )
  1236.     {
  1237.     if ( col < ccolso2 || col >= cols - ccolso2 )
  1238.         outputrow[col] = rowptr[crowso2][col];
  1239.     else if ( col == ccolso2 )
  1240.         {
  1241.         leftcol = col - ccolso2;
  1242.         for ( crow = 0; crow < crows; ++crow )
  1243.         {
  1244.         temprptr = rowptr[crow] + leftcol;
  1245.         for ( ccol = 0; ccol < ccols; ++ccol )
  1246.             {
  1247.             rcolumnsum[leftcol + ccol] += 
  1248.             PPM_GETR( *(temprptr + ccol) );
  1249.             gcolumnsum[leftcol + ccol] += 
  1250.             PPM_GETG( *(temprptr + ccol) );
  1251.             bcolumnsum[leftcol + ccol] += 
  1252.             PPM_GETB( *(temprptr + ccol) );
  1253.             }
  1254.         }
  1255.         for ( ccol = 0; ccol < ccols; ++ccol)
  1256.         {
  1257.         risum += rcolumnsum[leftcol + ccol];
  1258.         gisum += gcolumnsum[leftcol + ccol];
  1259.         bisum += bcolumnsum[leftcol + ccol];
  1260.         }
  1261.         temprsum = (float) risum * rmeanweight + 0.5;
  1262.         tempgsum = (float) gisum * gmeanweight + 0.5;
  1263.         tempbsum = (float) bisum * bmeanweight + 0.5;
  1264.         CHECK_RED;
  1265.         CHECK_GREEN;
  1266.         CHECK_BLUE;
  1267.         PPM_ASSIGN( outputrow[col], r, g, b );
  1268.         }
  1269.     else
  1270.         {
  1271.         /* Column numbers to subtract or add to isum */
  1272.         subcol = col - ccolso2 - 1;
  1273.         addcol = col + ccolso2;  
  1274.         for ( crow = 0; crow < crows; ++crow )
  1275.         {
  1276.         rcolumnsum[addcol] += PPM_GETR( rowptr[crow][addcol] );
  1277.         gcolumnsum[addcol] += PPM_GETG( rowptr[crow][addcol] );
  1278.         bcolumnsum[addcol] += PPM_GETB( rowptr[crow][addcol] );
  1279.         }
  1280.         risum = risum - rcolumnsum[subcol] + rcolumnsum[addcol];
  1281.         gisum = gisum - gcolumnsum[subcol] + gcolumnsum[addcol];
  1282.         bisum = bisum - bcolumnsum[subcol] + bcolumnsum[addcol];
  1283.         temprsum = (float) risum * rmeanweight + 0.5;
  1284.         tempgsum = (float) gisum * gmeanweight + 0.5;
  1285.         tempbsum = (float) bisum * bmeanweight + 0.5;
  1286.         CHECK_RED;
  1287.         CHECK_GREEN;
  1288.         CHECK_BLUE;
  1289.         PPM_ASSIGN( outputrow[col], r, g, b );
  1290.         }
  1291.     }
  1292.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  1293.  
  1294.     ++row;
  1295.     /* For all subsequent rows do it this way as the columnsums have been
  1296.     ** generated.  Now we can use them to reduce further calculations.
  1297.     */
  1298.     crowsp1 = crows + 1;
  1299.     for ( ; row < rows; ++row )
  1300.     {
  1301.     toprow = row + 1;
  1302.     temprow = row % (crows + 1);
  1303.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  1304.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1305.         pnm_promoteformatrow(
  1306.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  1307.  
  1308.     /* This rearrangement using crows+1 rowptrs and xelbufs will cause
  1309.     ** rowptr[0..crows-1] to always hold active xelbufs and for 
  1310.     ** rowptr[crows] to always hold the oldest (top most) xelbuf.
  1311.     */
  1312.     temprow = (toprow + 1) % crowsp1;
  1313.     i = 0;
  1314.     for (irow = temprow; irow < crowsp1; ++i, ++irow)
  1315.         rowptr[i] = xelbuf[irow];
  1316.     for (irow = 0; irow < temprow; ++irow, ++i)
  1317.         rowptr[i] = xelbuf[irow];
  1318.  
  1319.     risum = 0L;
  1320.     gisum = 0L;
  1321.     bisum = 0L;
  1322.     for ( col = 0; col < cols; ++col )
  1323.         {
  1324.         if ( col < ccolso2 || col >= cols - ccolso2 )
  1325.         outputrow[col] = rowptr[crowso2][col];
  1326.         else if ( col == ccolso2 )
  1327.         {
  1328.         leftcol = col - ccolso2;
  1329.         for ( ccol = 0; ccol < ccols; ++ccol )
  1330.             {
  1331.             tempcol = leftcol + ccol;
  1332.             rcolumnsum[tempcol] = rcolumnsum[tempcol]
  1333.             - PPM_GETR( rowptr[subrow][ccol] )
  1334.             + PPM_GETR( rowptr[addrow][ccol] );
  1335.             risum += rcolumnsum[tempcol];
  1336.             gcolumnsum[tempcol] = gcolumnsum[tempcol]
  1337.             - PPM_GETG( rowptr[subrow][ccol] )
  1338.             + PPM_GETG( rowptr[addrow][ccol] );
  1339.             gisum += gcolumnsum[tempcol];
  1340.             bcolumnsum[tempcol] = bcolumnsum[tempcol]
  1341.             - PPM_GETB( rowptr[subrow][ccol] )
  1342.             + PPM_GETB( rowptr[addrow][ccol] );
  1343.             bisum += bcolumnsum[tempcol];
  1344.             }
  1345.         temprsum = (float) risum * rmeanweight + 0.5;
  1346.         tempgsum = (float) gisum * gmeanweight + 0.5;
  1347.         tempbsum = (float) bisum * bmeanweight + 0.5;
  1348.         CHECK_RED;
  1349.         CHECK_GREEN;
  1350.         CHECK_BLUE;
  1351.         PPM_ASSIGN( outputrow[col], r, g, b );
  1352.         }
  1353.         else
  1354.         {
  1355.         /* Column numbers to subtract or add to isum */
  1356.         subcol = col - ccolso2 - 1;
  1357.         addcol = col + ccolso2;  
  1358.         rcolumnsum[addcol] = rcolumnsum[addcol]
  1359.             - PPM_GETR( rowptr[subrow][addcol] )
  1360.             + PPM_GETR( rowptr[addrow][addcol] );
  1361.         risum = risum - rcolumnsum[subcol] + rcolumnsum[addcol];
  1362.         gcolumnsum[addcol] = gcolumnsum[addcol]
  1363.             - PPM_GETG( rowptr[subrow][addcol] )
  1364.             + PPM_GETG( rowptr[addrow][addcol] );
  1365.         gisum = gisum - gcolumnsum[subcol] + gcolumnsum[addcol];
  1366.         bcolumnsum[addcol] = bcolumnsum[addcol]
  1367.             - PPM_GETB( rowptr[subrow][addcol] )
  1368.             + PPM_GETB( rowptr[addrow][addcol] );
  1369.         bisum = bisum - bcolumnsum[subcol] + bcolumnsum[addcol];
  1370.         temprsum = (float) risum * rmeanweight + 0.5;
  1371.         tempgsum = (float) gisum * gmeanweight + 0.5;
  1372.         tempbsum = (float) bisum * bmeanweight + 0.5;
  1373.         CHECK_RED;
  1374.         CHECK_GREEN;
  1375.         CHECK_BLUE;
  1376.         PPM_ASSIGN( outputrow[col], r, g, b );
  1377.         }
  1378.         }
  1379.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  1380.     }
  1381.  
  1382.     /* Now write out the remaining unconvolved rows in xelbuf. */
  1383.     for ( irow = crowso2 + 1; irow < crows; ++irow )
  1384.     pnm_writepnmrow(
  1385.             stdout, rowptr[irow], cols, maxval, newformat, 0 );
  1386.  
  1387.     pm_close( stdout );
  1388.     }
  1389.  
  1390.  
  1391. /* PPM Horizontal Convolution
  1392. **
  1393. ** Same as pgm_horizontal_convolve()
  1394. **
  1395. **/
  1396.  
  1397. static void
  1398. ppm_horizontal_convolve()
  1399.     {
  1400.     register int ccol, col;
  1401.     xel** xelbuf;
  1402.     xel* outputrow;
  1403.     xel x;
  1404.     xelval r, g, b;
  1405.     int row, crow;
  1406.     xel **rowptr, *temprptr;
  1407.     int leftcol;
  1408.     int i, irow;
  1409.     int temprow;
  1410.     int subcol, addcol;
  1411.     float rsum, gsum, bsum;
  1412.     int addrow, subrow;
  1413.     long **rrowsum, **rrowsumptr;
  1414.     long **growsum, **growsumptr;
  1415.     long **browsum, **browsumptr;
  1416.     int crowsp1;
  1417.     long temprsum, tempgsum, tempbsum;
  1418.  
  1419.     /* Allocate space for one convolution-matrix's worth of rows, plus
  1420.     ** a row output buffer. */
  1421.     xelbuf = pnm_allocarray( cols, crows + 1 );
  1422.     outputrow = pnm_allocrow( cols );
  1423.  
  1424.     /* Allocate array of pointers to xelbuf */
  1425.     rowptr = (xel **) pnm_allocarray( 1, crows + 1);
  1426.  
  1427.     /* Allocate intermediate row sums.  HORIZONTAL uses an extra row */
  1428.     rrowsum = (long **) pm_allocarray( cols, crows + 1, sizeof(float) );
  1429.     rrowsumptr = (long **) pnm_allocarray( 1, crows + 1);
  1430.     growsum = (long **) pm_allocarray( cols, crows + 1, sizeof(float) );
  1431.     growsumptr = (long **) pnm_allocarray( 1, crows + 1);
  1432.     browsum = (long **) pm_allocarray( cols, crows + 1, sizeof(float) );
  1433.     browsumptr = (long **) pnm_allocarray( 1, crows + 1);
  1434.  
  1435.     pnm_writepnminit( stdout, cols, rows, maxval, newformat, 0 );
  1436.  
  1437.     /* Read in one convolution-matrix's worth of image, less one row. */
  1438.     for ( row = 0; row < crows - 1; ++row )
  1439.     {
  1440.     pnm_readpnmrow( ifp, xelbuf[row], cols, maxval, format );
  1441.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1442.         pnm_promoteformatrow(
  1443.         xelbuf[row], cols, maxval, format, maxval, newformat );
  1444.     /* Write out just the part we're not going to convolve. */
  1445.     if ( row < crowso2 )
  1446.         pnm_writepnmrow( stdout, xelbuf[row], cols, maxval, newformat, 0 );
  1447.     }
  1448.  
  1449.     /* First row only */
  1450.     temprow = row % crows;
  1451.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  1452.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1453.     pnm_promoteformatrow(
  1454.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  1455.  
  1456.     temprow = (row + 1) % crows;
  1457.     i = 0;
  1458.     for (irow = temprow; irow < crows; ++i, ++irow)
  1459.     rowptr[i] = xelbuf[irow];
  1460.     for (irow = 0; irow < temprow; ++irow, ++i)
  1461.     rowptr[i] = xelbuf[irow];
  1462.  
  1463.     for ( crow = 0; crow < crows; ++crow )
  1464.     {
  1465.     rrowsumptr[crow] = rrowsum[crow];
  1466.     growsumptr[crow] = growsum[crow];
  1467.     browsumptr[crow] = browsum[crow];
  1468.     }
  1469.  
  1470.     for ( col = 0; col < cols; ++col )
  1471.     {
  1472.     if ( col < ccolso2 || col >= cols - ccolso2 )
  1473.         outputrow[col] = rowptr[crowso2][col];
  1474.     else if ( col == ccolso2 )
  1475.         {
  1476.         leftcol = col - ccolso2;
  1477.         rsum = 0.0;
  1478.         gsum = 0.0;
  1479.         bsum = 0.0;
  1480.         for ( crow = 0; crow < crows; ++crow )
  1481.         {
  1482.         temprptr = rowptr[crow] + leftcol;
  1483.         rrowsumptr[crow][leftcol] = 0L;
  1484.         growsumptr[crow][leftcol] = 0L;
  1485.         browsumptr[crow][leftcol] = 0L;
  1486.         for ( ccol = 0; ccol < ccols; ++ccol )
  1487.             {
  1488.             rrowsumptr[crow][leftcol] += 
  1489.                 PPM_GETR( *(temprptr + ccol) );
  1490.             growsumptr[crow][leftcol] += 
  1491.                 PPM_GETG( *(temprptr + ccol) );
  1492.             browsumptr[crow][leftcol] += 
  1493.                 PPM_GETB( *(temprptr + ccol) );
  1494.             }
  1495.         rsum += rrowsumptr[crow][leftcol] * rweights[crow][0];
  1496.         gsum += growsumptr[crow][leftcol] * gweights[crow][0];
  1497.         bsum += browsumptr[crow][leftcol] * bweights[crow][0];
  1498.         }
  1499.         temprsum = rsum + 0.5;
  1500.         tempgsum = gsum + 0.5;
  1501.         tempbsum = bsum + 0.5;
  1502.         CHECK_RED;
  1503.         CHECK_GREEN;
  1504.         CHECK_BLUE;
  1505.         PPM_ASSIGN( outputrow[col], r, g, b );
  1506.         }
  1507.     else
  1508.         {
  1509.         rsum = 0.0;
  1510.         gsum = 0.0;
  1511.         bsum = 0.0;
  1512.         leftcol = col - ccolso2;
  1513.         subcol = col - ccolso2 - 1;
  1514.         addcol = col + ccolso2;
  1515.         for ( crow = 0; crow < crows; ++crow )
  1516.         {
  1517.         rrowsumptr[crow][leftcol] = rrowsumptr[crow][subcol]
  1518.             - PPM_GETR( rowptr[crow][subcol] )
  1519.             + PPM_GETR( rowptr[crow][addcol] );
  1520.         rsum += rrowsumptr[crow][leftcol] * rweights[crow][0];
  1521.         growsumptr[crow][leftcol] = growsumptr[crow][subcol]
  1522.             - PPM_GETG( rowptr[crow][subcol] )
  1523.             + PPM_GETG( rowptr[crow][addcol] );
  1524.         gsum += growsumptr[crow][leftcol] * gweights[crow][0];
  1525.         browsumptr[crow][leftcol] = browsumptr[crow][subcol]
  1526.             - PPM_GETB( rowptr[crow][subcol] )
  1527.             + PPM_GETB( rowptr[crow][addcol] );
  1528.         bsum += browsumptr[crow][leftcol] * bweights[crow][0];
  1529.         }
  1530.         temprsum = rsum + 0.5;
  1531.         tempgsum = gsum + 0.5;
  1532.         tempbsum = bsum + 0.5;
  1533.         CHECK_RED;
  1534.         CHECK_GREEN;
  1535.         CHECK_BLUE;
  1536.         PPM_ASSIGN( outputrow[col], r, g, b );
  1537.         }
  1538.         }
  1539.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  1540.  
  1541.  
  1542.     /* For all subsequent rows */
  1543.  
  1544.     subrow = crows;
  1545.     addrow = crows - 1;
  1546.     crowsp1 = crows + 1;
  1547.     ++row;
  1548.     for ( ; row < rows; ++row )
  1549.     {
  1550.     temprow = row % crowsp1;
  1551.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  1552.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1553.         pnm_promoteformatrow(
  1554.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  1555.  
  1556.     temprow = (row + 2) % crowsp1;
  1557.     i = 0;
  1558.     for (irow = temprow; irow < crowsp1; ++i, ++irow)
  1559.         {
  1560.         rowptr[i] = xelbuf[irow];
  1561.         rrowsumptr[i] = rrowsum[irow];
  1562.         growsumptr[i] = growsum[irow];
  1563.         browsumptr[i] = browsum[irow];
  1564.         }
  1565.     for (irow = 0; irow < temprow; ++irow, ++i)
  1566.         {
  1567.         rowptr[i] = xelbuf[irow];
  1568.         rrowsumptr[i] = rrowsum[irow];
  1569.         growsumptr[i] = growsum[irow];
  1570.         browsumptr[i] = browsum[irow];
  1571.         }
  1572.  
  1573.     for ( col = 0; col < cols; ++col )
  1574.         {
  1575.         if ( col < ccolso2 || col >= cols - ccolso2 )
  1576.         outputrow[col] = rowptr[crowso2][col];
  1577.         else if ( col == ccolso2 )
  1578.         {
  1579.         rsum = 0.0;
  1580.         gsum = 0.0;
  1581.         bsum = 0.0;
  1582.         leftcol = col - ccolso2;
  1583.         rrowsumptr[addrow][leftcol] = 0L;
  1584.         growsumptr[addrow][leftcol] = 0L;
  1585.         browsumptr[addrow][leftcol] = 0L;
  1586.         for ( ccol = 0; ccol < ccols; ++ccol )
  1587.             {
  1588.             rrowsumptr[addrow][leftcol] += 
  1589.             PPM_GETR( rowptr[addrow][leftcol + ccol] );
  1590.             growsumptr[addrow][leftcol] += 
  1591.             PPM_GETG( rowptr[addrow][leftcol + ccol] );
  1592.             browsumptr[addrow][leftcol] += 
  1593.             PPM_GETB( rowptr[addrow][leftcol + ccol] );
  1594.             }
  1595.         for ( crow = 0; crow < crows; ++crow )
  1596.             {
  1597.             rsum += rrowsumptr[crow][leftcol] * rweights[crow][0];
  1598.             gsum += growsumptr[crow][leftcol] * gweights[crow][0];
  1599.             bsum += browsumptr[crow][leftcol] * bweights[crow][0];
  1600.             }
  1601.         temprsum = rsum + 0.5;
  1602.         tempgsum = gsum + 0.5;
  1603.         tempbsum = bsum + 0.5;
  1604.         CHECK_RED;
  1605.         CHECK_GREEN;
  1606.         CHECK_BLUE;
  1607.         PPM_ASSIGN( outputrow[col], r, g, b );
  1608.         }
  1609.         else
  1610.         {
  1611.         rsum = 0.0;
  1612.         gsum = 0.0;
  1613.         bsum = 0.0;
  1614.         leftcol = col - ccolso2;
  1615.         subcol = col - ccolso2 - 1;
  1616.         addcol = col + ccolso2;  
  1617.         rrowsumptr[addrow][leftcol] = rrowsumptr[addrow][subcol]
  1618.             - PPM_GETR( rowptr[addrow][subcol] )
  1619.             + PPM_GETR( rowptr[addrow][addcol] );
  1620.         growsumptr[addrow][leftcol] = growsumptr[addrow][subcol]
  1621.             - PPM_GETG( rowptr[addrow][subcol] )
  1622.             + PPM_GETG( rowptr[addrow][addcol] );
  1623.         browsumptr[addrow][leftcol] = browsumptr[addrow][subcol]
  1624.             - PPM_GETB( rowptr[addrow][subcol] )
  1625.             + PPM_GETB( rowptr[addrow][addcol] );
  1626.         for ( crow = 0; crow < crows; ++crow )
  1627.             {
  1628.             rsum += rrowsumptr[crow][leftcol] * rweights[crow][0];
  1629.             gsum += growsumptr[crow][leftcol] * gweights[crow][0];
  1630.             bsum += browsumptr[crow][leftcol] * bweights[crow][0];
  1631.             }
  1632.         temprsum = rsum + 0.5;
  1633.         tempgsum = gsum + 0.5;
  1634.         tempbsum = bsum + 0.5;
  1635.         CHECK_RED;
  1636.         CHECK_GREEN;
  1637.         CHECK_BLUE;
  1638.         PPM_ASSIGN( outputrow[col], r, g, b );
  1639.         }
  1640.         }
  1641.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  1642.     }
  1643.  
  1644.     /* Now write out the remaining unconvolved rows in xelbuf. */
  1645.     for ( irow = crowso2 + 1; irow < crows; ++irow )
  1646.     pnm_writepnmrow(
  1647.             stdout, rowptr[irow], cols, maxval, newformat, 0 );
  1648.  
  1649.     pm_close( stdout );
  1650.     }
  1651.  
  1652.  
  1653. /* PPM Vertical Convolution
  1654. **
  1655. ** Same as pgm_vertical_convolve()
  1656. **
  1657. */
  1658.  
  1659. static void
  1660. ppm_vertical_convolve()
  1661.     {
  1662.     register int ccol, col;
  1663.     xel** xelbuf;
  1664.     xel* outputrow;
  1665.     xelval r, g, b;
  1666.     int row, crow;
  1667.     xel **rowptr, *temprptr;
  1668.     int leftcol;
  1669.     int i, irow;
  1670.     int toprow, temprow;
  1671.     int subrow, addrow;
  1672.     int tempcol;
  1673.     float rsum, gsum, bsum;
  1674.     long *rcolumnsum, *gcolumnsum, *bcolumnsum;
  1675.     int crowsp1;
  1676.     int addcol;
  1677.     long temprsum, tempgsum, tempbsum;
  1678.  
  1679.     /* Allocate space for one convolution-matrix's worth of rows, plus
  1680.     ** a row output buffer. VERTICAL uses an extra row. */
  1681.     xelbuf = pnm_allocarray( cols, crows + 1 );
  1682.     outputrow = pnm_allocrow( cols );
  1683.  
  1684.     /* Allocate array of pointers to xelbuf */
  1685.     rowptr = (xel **) pnm_allocarray( 1, crows + 1 );
  1686.  
  1687.     /* Allocate space for intermediate column sums */
  1688.     rcolumnsum = (long *) pm_allocrow( cols, sizeof(long) );
  1689.     gcolumnsum = (long *) pm_allocrow( cols, sizeof(long) );
  1690.     bcolumnsum = (long *) pm_allocrow( cols, sizeof(long) );
  1691.     for ( col = 0; col < cols; ++col )
  1692.     {
  1693.     rcolumnsum[col] = 0L;
  1694.     gcolumnsum[col] = 0L;
  1695.     bcolumnsum[col] = 0L;
  1696.     }
  1697.  
  1698.     pnm_writepnminit( stdout, cols, rows, maxval, newformat, 0 );
  1699.  
  1700.     /* Read in one convolution-matrix's worth of image, less one row. */
  1701.     for ( row = 0; row < crows - 1; ++row )
  1702.     {
  1703.     pnm_readpnmrow( ifp, xelbuf[row], cols, maxval, format );
  1704.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1705.         pnm_promoteformatrow(
  1706.         xelbuf[row], cols, maxval, format, maxval, newformat );
  1707.     /* Write out just the part we're not going to convolve. */
  1708.     if ( row < crowso2 )
  1709.         pnm_writepnmrow( stdout, xelbuf[row], cols, maxval, newformat, 0 );
  1710.     }
  1711.  
  1712.     /* Now the rest of the image - read in the row at the end of
  1713.     ** xelbuf, and convolve and write out the row in the middle.
  1714.     */
  1715.     /* For first row only */
  1716.  
  1717.     toprow = row + 1;
  1718.     temprow = row % crows;
  1719.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  1720.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1721.     pnm_promoteformatrow(
  1722.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  1723.  
  1724.     /* Arrange rowptr to eliminate the use of mod function to determine
  1725.     ** which row of xelbuf is 0...crows.  Mod function can be very costly.
  1726.     */
  1727.     temprow = toprow % crows;
  1728.     i = 0;
  1729.     for (irow = temprow; irow < crows; ++i, ++irow)
  1730.     rowptr[i] = xelbuf[irow];
  1731.     for (irow = 0; irow < temprow; ++irow, ++i)
  1732.     rowptr[i] = xelbuf[irow];
  1733.  
  1734.     for ( col = 0; col < cols; ++col )
  1735.     {
  1736.     if ( col < ccolso2 || col >= cols - ccolso2 )
  1737.         outputrow[col] = rowptr[crowso2][col];
  1738.     else if ( col == ccolso2 )
  1739.         {
  1740.         rsum = 0.0;
  1741.         gsum = 0.0;
  1742.         bsum = 0.0;
  1743.         leftcol = col - ccolso2;
  1744.         for ( crow = 0; crow < crows; ++crow )
  1745.         {
  1746.         temprptr = rowptr[crow] + leftcol;
  1747.         for ( ccol = 0; ccol < ccols; ++ccol )
  1748.             {
  1749.             rcolumnsum[leftcol + ccol] += 
  1750.             PPM_GETR( *(temprptr + ccol) );
  1751.             gcolumnsum[leftcol + ccol] += 
  1752.             PPM_GETG( *(temprptr + ccol) );
  1753.             bcolumnsum[leftcol + ccol] += 
  1754.             PPM_GETB( *(temprptr + ccol) );
  1755.             }
  1756.         }
  1757.         for ( ccol = 0; ccol < ccols; ++ccol)
  1758.         {
  1759.         rsum += rcolumnsum[leftcol + ccol] * rweights[0][ccol];
  1760.         gsum += gcolumnsum[leftcol + ccol] * gweights[0][ccol];
  1761.         bsum += bcolumnsum[leftcol + ccol] * bweights[0][ccol];
  1762.         }
  1763.         temprsum = rsum + 0.5;
  1764.         tempgsum = gsum + 0.5;
  1765.         tempbsum = bsum + 0.5;
  1766.         CHECK_RED;
  1767.         CHECK_GREEN;
  1768.         CHECK_BLUE;
  1769.         PPM_ASSIGN( outputrow[col], r, g, b );
  1770.         }
  1771.     else
  1772.         {
  1773.         rsum = 0.0;
  1774.         gsum = 0.0;
  1775.         bsum = 0.0;
  1776.         leftcol = col - ccolso2;
  1777.         addcol = col + ccolso2;  
  1778.         for ( crow = 0; crow < crows; ++crow )
  1779.         {
  1780.         rcolumnsum[addcol] += PPM_GETR( rowptr[crow][addcol] );
  1781.         gcolumnsum[addcol] += PPM_GETG( rowptr[crow][addcol] );
  1782.         bcolumnsum[addcol] += PPM_GETB( rowptr[crow][addcol] );
  1783.         }
  1784.         for ( ccol = 0; ccol < ccols; ++ccol )
  1785.         {
  1786.         rsum += rcolumnsum[leftcol + ccol] * rweights[0][ccol];
  1787.         gsum += gcolumnsum[leftcol + ccol] * gweights[0][ccol];
  1788.         bsum += bcolumnsum[leftcol + ccol] * bweights[0][ccol];
  1789.         }
  1790.         temprsum = rsum + 0.5;
  1791.         tempgsum = gsum + 0.5;
  1792.         tempbsum = bsum + 0.5;
  1793.         CHECK_RED;
  1794.         CHECK_GREEN;
  1795.         CHECK_BLUE;
  1796.         PPM_ASSIGN( outputrow[col], r, g, b );
  1797.         }
  1798.     }
  1799.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  1800.  
  1801.     /* For all subsequent rows */
  1802.     subrow = crows;
  1803.     addrow = crows - 1;
  1804.     crowsp1 = crows + 1;
  1805.     ++row;
  1806.     for ( ; row < rows; ++row )
  1807.     {
  1808.     toprow = row + 1;
  1809.     temprow = row % (crows +1);
  1810.     pnm_readpnmrow( ifp, xelbuf[temprow], cols, maxval, format );
  1811.     if ( PNM_FORMAT_TYPE(format) != newformat )
  1812.         pnm_promoteformatrow(
  1813.         xelbuf[temprow], cols, maxval, format, maxval, newformat );
  1814.  
  1815.     /* Arrange rowptr to eliminate the use of mod function to determine
  1816.     ** which row of xelbuf is 0...crows.  Mod function can be very costly.
  1817.     */
  1818.     temprow = (toprow + 1) % crowsp1;
  1819.     i = 0;
  1820.     for (irow = temprow; irow < crowsp1; ++i, ++irow)
  1821.         rowptr[i] = xelbuf[irow];
  1822.     for (irow = 0; irow < temprow; ++irow, ++i)
  1823.         rowptr[i] = xelbuf[irow];
  1824.  
  1825.     for ( col = 0; col < cols; ++col )
  1826.         {
  1827.         if ( col < ccolso2 || col >= cols - ccolso2 )
  1828.         outputrow[col] = rowptr[crowso2][col];
  1829.         else if ( col == ccolso2 )
  1830.         {
  1831.         rsum = 0.0;
  1832.         gsum = 0.0;
  1833.         bsum = 0.0;
  1834.         leftcol = col - ccolso2;
  1835.         for ( ccol = 0; ccol < ccols; ++ccol )
  1836.             {
  1837.             tempcol = leftcol + ccol;
  1838.             rcolumnsum[tempcol] = rcolumnsum[tempcol] 
  1839.             - PPM_GETR( rowptr[subrow][ccol] )
  1840.             + PPM_GETR( rowptr[addrow][ccol] );
  1841.             rsum = rsum + rcolumnsum[tempcol] * rweights[0][ccol];
  1842.             gcolumnsum[tempcol] = gcolumnsum[tempcol] 
  1843.             - PPM_GETG( rowptr[subrow][ccol] )
  1844.             + PPM_GETG( rowptr[addrow][ccol] );
  1845.             gsum = gsum + gcolumnsum[tempcol] * gweights[0][ccol];
  1846.             bcolumnsum[tempcol] = bcolumnsum[tempcol] 
  1847.             - PPM_GETB( rowptr[subrow][ccol] )
  1848.             + PPM_GETB( rowptr[addrow][ccol] );
  1849.             bsum = bsum + bcolumnsum[tempcol] * bweights[0][ccol];
  1850.             }
  1851.         temprsum = rsum + 0.5;
  1852.         tempgsum = gsum + 0.5;
  1853.         tempbsum = bsum + 0.5;
  1854.         CHECK_RED;
  1855.         CHECK_GREEN;
  1856.         CHECK_BLUE;
  1857.         PPM_ASSIGN( outputrow[col], r, g, b );
  1858.         }
  1859.         else
  1860.         {
  1861.         rsum = 0.0;
  1862.         gsum = 0.0;
  1863.         bsum = 0.0;
  1864.         leftcol = col - ccolso2;
  1865.         addcol = col + ccolso2;
  1866.         rcolumnsum[addcol] = rcolumnsum[addcol]
  1867.             - PPM_GETR( rowptr[subrow][addcol] )
  1868.             + PPM_GETR( rowptr[addrow][addcol] );
  1869.         gcolumnsum[addcol] = gcolumnsum[addcol]
  1870.             - PPM_GETG( rowptr[subrow][addcol] )
  1871.             + PPM_GETG( rowptr[addrow][addcol] );
  1872.         bcolumnsum[addcol] = bcolumnsum[addcol]
  1873.             - PPM_GETB( rowptr[subrow][addcol] )
  1874.             + PPM_GETB( rowptr[addrow][addcol] );
  1875.         for ( ccol = 0; ccol < ccols; ++ccol )
  1876.             {
  1877.             rsum += rcolumnsum[leftcol + ccol] * rweights[0][ccol];
  1878.             gsum += gcolumnsum[leftcol + ccol] * gweights[0][ccol];
  1879.             bsum += bcolumnsum[leftcol + ccol] * bweights[0][ccol];
  1880.             }
  1881.         temprsum = rsum + 0.5;
  1882.         tempgsum = gsum + 0.5;
  1883.         tempbsum = bsum + 0.5;
  1884.         CHECK_RED;
  1885.         CHECK_GREEN;
  1886.         CHECK_BLUE;
  1887.         PPM_ASSIGN( outputrow[col], r, g, b );
  1888.         }
  1889.         }
  1890.     pnm_writepnmrow( stdout, outputrow, cols, maxval, newformat, 0 );
  1891.     }
  1892.  
  1893.     /* Now write out the remaining unconvolved rows in xelbuf. */
  1894.     for ( irow = crowso2 + 1; irow < crows; ++irow )
  1895.     pnm_writepnmrow(
  1896.             stdout, rowptr[irow], cols, maxval, newformat, 0 );
  1897.  
  1898.     pm_close( stdout );
  1899.     }
  1900.