C
C     Open Systems Laboratory
C     http://www.lam-mpi.org/tutorials/
C     Indiana University
C
C     MPI Tutorial
C     Lab : Image processing -- sum of squares
C
C     Mail questions regarding tutorial material to lam at lam dash mpi dot org
C
      program main
      include 'mpif.h'

      integer   iwidth, iheight, numpixels
      integer   rank, comm_size, sum, my_sum
      integer   i, my_count, ierr
      real      rms
      
      character recvbuf(65536)
      character pixels(65536)

      call MPI_INIT(ierr)
      call MPI_COMM_RANK(MPI_COMM_WORLD, rank, ierr)
      call MPI_COMM_SIZE(MPI_COMM_WORLD, comm_size, ierr)

      if (rank .eq. 0) then

C        Load the image

         iheight = 500
         iwidth = 500
         call LOADTIFF('irish.tif', pixels, iwidth, iheight)
         numpixels = iwidth * iheight

         my_count = numpixels / comm_size
      endif

C     Broadcasts my_count to all the processes 

      call MPI_BCAST(my_count, 1, MPI_INTEGER, 0, 
     $     MPI_COMM_WORLD, ierr)

C     Scatter the image 

      call MPI_SCATTER(pixels, my_count, MPI_CHARACTER, recvbuf, 
     $     my_count, MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr)

C     Take the sum of the squares of the partial image

      my_sum = 0
      do i = 1, my_count
         my_sum = my_sum + ichar(recvbuf(i)) * ichar(recvbuf(i))
      enddo

C     Find the global sum of the squares

      call MPI_REDUCE(my_sum, sum, 1, MPI_INTEGER, MPI_SUM, 0, 
     $     MPI_COMM_WORLD, ierr)
      
C     Rank 0 calculates the root mean square 

      if (rank .eq. 0) then
         rms = sqrt(real(sum) / real(numpixels))
         print *, 'RMS = ', rms
      endif
      
      call MPI_FINALIZE(ierr)
      stop
      end