RECURSIVE SUBROUTINE findParticlesInCell(cellNumber, bnd_box_min, bnd_box_max, blockCells)

	USE fly_h

        INCLUDE 'mpif.h'
	integer(kind=4) :: cellNumber
	real, dimension(ndim) :: bnd_box_min, bnd_box_max
	real(KIND=8), dimension(nZones(1), nZones(2), nZones(3)) :: blockCells

#ifdef FLASH

	integer :: k,ind_pe
	integer(kind=4), dimension(nsubcell) :: isub
	integer, dimension(ndim) :: cellIndex

	INTEGER(KIND=MPI_ADDRESS_KIND) :: startIndex

	real :: blockSize
	real, dimension(ndim) :: cellSize
	real(KIND=8), dimension(ndim) :: bodyPosition

!===================================================================

        i_sh = cellNumber - nbodsmax

	
	IF(POW2) THEN
        ind_pe  = IAND((i_sh - 1), npes - 1)
        ind_loc = ISHFT((i_sh - 1), -lpes) + 1
	ELSE
        ind_pe  = MOD((i_sh - 1), npes)
        ind_loc =((i_sh - 1)/ npes) + 1
	ENDIF
! 	write(100+me,*)"ME=",me," cella_number=",cellNumber," ind_pe=",ind_pe
!	call flush(uterm)

        IF(ind_pe.EQ.me) THEN

	   isub(1:nsubcell)	= subp(1:8, ind_loc)

        ELSE
! 	write(100+me,*)"ME=",me," passo1 isub=",isub
!	call flush(uterm)

	   startIndex = (((ind_loc - 1) * nsubcell) + 1) - 1

           CALL MPI_WIN_LOCK(MPI_LOCK_SHARED, ind_pe, 0, win_subp, ierror)
           CALL MPI_GET(isub(1), nsubcell, MPI_INTEGER4, ind_pe,        &
     &       		startIndex, nsubcell, MPI_INTEGER4, win_subp, ierror)
           CALL MPI_WIN_UNLOCK(ind_pe, win_subp, ierror)
!	write(uterm,*)"ME=",me," passo2"
!	call flush(uterm)

        ENDIF
!	write(100+me,*)"ME=",me," isub=",isub
!	call flush(uterm)

	blockSize       = bnd_box_max(1) - bnd_box_min(1)
        cellSize	= blockSize / nZones

	DO k = 1, nsubcell  ! QUI controllare che arrivi dove dovrebbe arrivare.
	
	   IF(isub(k).GT.nbodsmax) THEN

!	write(uterm,*)"ME=",me," passo3"
!	call flush(uterm)
	      CALL findParticlesInCell(isub(k), bnd_box_min, &
					bnd_box_max, blockCells) 

	   ELSEIF(isub(k).GT.0) THEN

	      i_sh	= isub(k)
	      ind_pe	= (i_sh-1)/ max_pr
	      ind_loc	= i_sh - (max_pr * ind_pe)

	      IF(ind_pe.EQ.me) THEN
			
	 	 bodyPosition(1:3) = pos(1:3, ind_loc)
				
	      ELSE

                 startIndex = (((ind_loc - 1) * ndim) + 1) - 1

!	write(uterm,*)"ME=",me," passo4"
!	call flush(uterm)
                 CALL MPI_WIN_LOCK(MPI_LOCK_SHARED, ind_pe, 0,          &
                                  win_pos, ierror)
                 CALL MPI_GET(bodyPosition(1), ndim, MPI_REAL8, ind_pe, &
                             startIndex, ndim, MPI_REAL8, win_pos, ierror)
                 CALL MPI_WIN_UNLOCK(ind_pe, win_pos, ierror) 

!	write(uterm,*)"ME=",me," passo5"
!	call flush(uterm)
	      ENDIF
!	IF(me.eq.1) write(100+me,*)"ME=",me," bodyPosition =",bodyPosition," bnd_box_min =",bnd_box_min, &
!	" bnd_box_max =",bnd_box_max
!	call flush(uterm)

	      IF(bodyPosition(3).GE.bnd_box_min(3).AND.bodyPosition(3).LT.bnd_box_max(3)) THEN
	      IF(bodyPosition(2).GE.bnd_box_min(2).AND.bodyPosition(2).LT.bnd_box_max(2)) THEN
	      IF(bodyPosition(1).GE.bnd_box_min(1).AND.bodyPosition(1).LT.bnd_box_max(1)) THEN

	      cellIndex = (bodyPosition - bnd_box_min)/cellSize + 1

              blockCells(cellIndex(1), cellIndex(2), cellIndex(3)) = &
				blockCells(cellIndex(1), cellIndex(2), cellIndex(3)) + 1
			
!	write(100+me,*)"ME=",me," OK cellIndex=",cellIndex
!	call flush(uterm)
	      ENDIF
	      ENDIF
	      ENDIF
	   ENDIF
	ENDDO
#endif
!	write(uterm,*)"ME=",me," esco da findParticlesInCell"
!	call flush(uterm)
	RETURN

	END