program mergebloubs

          !-------------------------------------------!
          !       THIS IS A BIG MESS OF BUGS          !
          !-------------------------------------------!

  use bloubspace
  implicit none

  integer, parameter       :: NB_MAX_BLOUBS = 25000

  character(200)           :: infile, outfile
  type(t_bloubs)           :: bloub, newbloub
  integer                  :: inu, outu, errcode

  type(t_bloubs), dimension(:), allocatable  :: les_bloubs
  integer                                    :: i, idx, nbr_merge
  real                                       :: rval
  logical                                    :: merged

! --------------- check command line parameters
  if (IARGC() .ne. 2) then
    STOP ": NEED IN AND OUT FILENAME"
  endif
  call getarg(1, infile)
  call getarg(2, outfile)

  write(0, '(A, 2A20, I8)') "*** mergebloubs ",        &
                 trim(infile), trim(outfile), NB_MAX_BLOUBS

! --------------- allocate memory for the people

  allocate (les_bloubs(NB_MAX_BLOUBS), stat=errcode)
  if (0 .NE. errcode) then
    STOP " : NO ENOUGH MEMORY"
  endif

  do i = 1, NB_MAX_BLOUBS
    les_bloubs(i)%alive = .FALSE.
  enddo

! --------------- open / creat the files

  open(newunit=inu,                                          &
                  file=trim(infile), form='unformatted',     &
                  iostat=errcode,                            &
                  action='read', status='old')
  if (0 .ne. errcode) then
    STOP " : CAN'T OPEN FILE " // trim(infile)
  endif

  open(newunit=outu,                                         &
                  file=trim(outfile), form='unformatted',    &
                  iostat=errcode,                            &
                  action='write', status='replace')
  if (0 .ne. errcode) then
    STOP " : CAN'T OPEN " // trim(outfile) // "FOR WRITE"
  endif

! --------------- read the first bloub

  idx = 1

  read (unit=inu, iostat=errcode) bloub
    if (0 .ne. errcode) then
      STOP " : ERR READING FIRST BLOUB"
    endif
  call display_bloub (bloub, "first bloub")
  write(outu, iostat=errcode) bloub
  if (0 .ne. errcode) then
    STOP " : FIRST BLOUB, WRITE ERROR TO " // trim(outfile)
  endif

  les_bloubs(idx) = bloub
  idx = idx + 1

! --------------- loop over the other bloubs
  nbr_merge = 0

  do                 ! infinite loop

    ! print *, "============ PASS ", idx
    if (idx .EQ. NB_MAX_BLOUBS) then
      write(0, '(I8, A)') idx, " max number of bloubs reached"
      exit
    endif

    ! read the next bloub from input file
    read (unit=inu, iostat=errcode) bloub
      if (0 .ne. errcode) then
        exit
      endif
    !! call display_bloub (bloub, "next bloub")

    if (.NOT. bloub%alive) then
      STOP " : I HAVE READ A DEAD BLOUB"
    endif

    ! check with all the previuous blobs
    merged = .FALSE.
    do i = 1, idx-1
      if (.NOT. les_bloubs(i)%alive) then
        ! print *, "dead bloub at ", i, " on ", idx
        ! call display_bloub(les_bloubs(i), "DEAD ? WTF ?")
        cycle
      endif
      rval = distance_of_bloubs(les_bloubs(i), bloub)
      if (rval .LT. (les_bloubs(i)%radius + bloub%radius)) then
        print *, "contact : ", i, idx, rval
        call merge_two_bloubs(les_bloubs(i), bloub, newbloub)
        les_bloubs(i)%alive = .FALSE.
        nbr_merge = nbr_merge + 1
        merged = .TRUE.
      endif
    enddo

    if (merged) then
      les_bloubs(idx) = newbloub
      bloub = newbloub
    else
      ! put old bloub in the list
      les_bloubs(idx) = bloub
    endif

    write(outu, iostat=errcode) bloub
    if (0 .ne. errcode) then
      STOP " : WRITE ERROR TO " // trim(outfile)
    endif

    idx = idx + 1

    ! print *, "idx = ", idx

  enddo              ! end of infinit... WHAT?

! --------------- is the job done ?

  close(inu)    ;    close(outu)

  write(0, '()')
  write(0, '(I5, A)') nbr_merge, " merges"
  write(0, '(A)') "--- mergebloubs . . . . . . .  [done]"

! --------------------------------------------------------------
contains

  subroutine merge_two_bloubs(bla, blb, blr)
      type(t_bloubs), intent(in)   :: bla, blb
      type(t_bloubs), intent(out)  :: blr

      blr%nick     = "newbie  "
      blr%num      = 0                !    ???

      blr%px = (bla%px + blb%px) / 2.0
      blr%py = (bla%py + blb%py) / 2.0
      blr%pz = (bla%pz + blb%pz) / 2.0

      blr%vx = (bla%vx + blb%vx) / 2.0
      blr%vy = (bla%vy + blb%vy) / 2.0
      blr%vz = (bla%vz + blb%vz) / 2.0

      blr%radius = (bla%radius + blb%radius) / 2.718
      blr%age = min(bla%age, blb%age)

      ! bring it to life !
      blr%alive    = .TRUE.

  end subroutine

end program