!***********************************************************************
!
!
                           SUBROUTINE read_s_bin
!
!
!***********************************************************************
!
!
!     Subroutine to read in the data associated with the bodies from
!     the BINARY  input file.  
!=======================================================================

	 USE fly_h 
         implicit none
	INCLUDE 'mpif.h'

!   Declaration of local variables.
!   -------------------------------
	INTEGER(KIND=MPI_ADDRESS_KIND) :: startIndex
        INTEGER  :: tag,istatus(MPI_STATUS_SIZE),i,pint
        INTEGER*8  :: p,icounttr,periodX,periodY,periodZ
        INTEGER*8 :: pbefore,body4Plane,residualBodyY,residualBodyX
        INTEGER, DIMENSION(N_PES*2) :: seed
        CHARACTER(LEN=4) :: ifl      
	CHARACTER(LEN=256) :: name_file
	INTEGER*8 :: NLONG8,itest
		REAL(KIND=8) :: CBRT,stepsize,stepsize4,stepsize8
		REAL :: X,Y,Z,RN,RN0
	REAL(KIND=8) ::c0a, c0b,  c1a, c1b, c2a,c2b
		 
!
! mass_read is the value of the mass of each body of the simulation
!
!=======================================================================
!   propagate the mass_read value to all bodies
	
	
	DO p=1,nb_loc
            mass(p)=mass_read
        ENDDO

!   Read in body data.
!   ------------------
        
	IF(nbodies.GT.nbodsmax)                                         &
          CALL error('read_b_bin error1: inconsistent input')
	

!   Compute numbers of intervals in the box
!   --------------------------
	CBRT=nbodsmax**(1.0/3.0)
	stepsize=rsize/CBRT
	stepsize4=stepsize/4.0
	stepsize8=stepsize/8.0
	DO i=1,N_PES*2
		seed(i)=me
	ENDDO
	CALL RANDOM_SEED(PUT=seed)
	pbefore=0
	DO PE_indx=1,me
		pbefore=pbefore+nb_res_loc(PE_indx)
	ENDDO
	IF (me.EQ.0) write(uterm,*) "PE=",me," Simulated bodies stepsize=",stepsize
!   --------------------------
	body4Plane=CBRT*CBRT
	periodZ=(pbefore/body4Plane)
	residualBodyY=pbefore-periodZ*body4Plane
	periodY=(residualBodyY/CBRT)
	residualBodyX=(residualBodyY-periodY*CBRT)
	periodX=residualBodyX
!   --------------------------
	   c0a=MPI_WTIME()
	  write(uterm,*) "TP0 PE=",me,CBRT,pbefore,body4Plane,periodZ,residualBodyY,periodY,residualBodyX,periodX
	  call flush(uterm)
!	  IF(me.gt.0) THEN
!	  	periodZ=periodZ-1
!	  	periodY=periodY-1
!	  	periodX=periodX-1
!	ENDIF
	  write(uterm,*) "TP0 PE=",me," Simulated bodies  Start", periodX,periodY,periodZ
	  call flush(uterm)

	   X=0
	   Y=0
	   Z=0
	   DO p=1,nb_res_loc(me+1)
        CALL RANDOM_NUMBER(RN0)
	    RN=RN0*stepsize4-stepsize8
	    X=(periodX)*stepsize+RN 
	    periodX=periodX+1
	    Y=(periodY)*stepsize+RN 
	    Z=(periodZ)*stepsize+RN 
	    IF(periodX.ge.CBRT) THEN
	    	periodX=0
	    	periodY=periodY+1
	    ENDIF
	    IF(periodY.ge.CBRT) THEN
	    	periodX=0
	    	periodY=0
	    	periodZ=periodZ+1
	    ENDIF
	    IF(periodZ.ge.CBRT) THEN
	    	periodX=0
	    	periodY=0
	    	periodZ=0
	    ENDIF
	    pos(1,p)=X
	    if(X.lt.0.0) pos(1,p)=0+stepsize8*RN0
	    if(X.gt.rsize) pos(1,p)=rsize-stepsize8*RN0
	    pos(2,p)=Y
	    if(Y.lt.0.0) pos(2,p)=0+stepsize8*RN0
	    if(Y.gt.rsize) pos(2,p)=rsize-stepsize8*RN0
	    pos(3,p)=Z
	    if(Z.lt.0.0) pos(3,p)=0+stepsize8*RN0
	    if(Z.gt.rsize) pos(3,p)=rsize-stepsize8*RN0
	  ENDDO 		

	  write(uterm,*) "TP0 PE=",me," Simulated bodies  END", periodX,periodY,periodZ
	  call flush(uterm)
	  
        
	
	
	CALL MPI_BARRIER(MPI_COMM_WORLD,ierror)

	  write(uterm,*) "PE=",me," Simulated vel for PE=",me," START"
	  call flush(uterm)
	   DO p=1,nb_res_loc(me+1)
	   		CALL RANDOM_NUMBER(RN0)
			vel(1,p)=RN0*200.0-100.0	
	   		CALL RANDOM_NUMBER(RN0)
			vel(2,p)=RN0*200.0-100.0	
	   		CALL RANDOM_NUMBER(RN0)
			vel(3,p)=RN0*200.0-100.0
	  ENDDO	
	  write(uterm,*) "PE=",me," Simulated vel for PE=",me," END"
	  call flush(uterm)
	  
        

	
	CALL MPI_BARRIER(MPI_COMM_WORLD,ierror)		

	IF(me.eq.0) THEN
	write(uterm,*) 'pos first = ', pos(1:3,1)
	write(uterm,*) 'vel first = ', vel(1:3,1)
	call flush(uterm)
	ENDIF
	
	CALL MPI_BARRIER(MPI_COMM_WORLD,ierror)

	IF(me.eq.npes-1) THEN
	write(uterm,*) 'pos last  = ', pos(1:3,nb_res_loc(me+1))
	write(uterm,*) 'vel last  = ', vel(1:3,nb_res_loc(me+1))
	call flush(uterm)
	ENDIF
	
	CALL MPI_BARRIER(MPI_COMM_WORLD,ierror)
	if(inasc.eq.'T') CALL write_b_asc
	if(inasc.eq.'R') CALL write_b_bin
	IF((inasc.eq.'R') .or. (inasc .eq. 'T')) THEN 
		if(me.eq.0) write(uterm,*) "Simulated output file created and stop"
		call MPI_Finalize(ierror)
		STOP 
	ENDIF
        
	RETURN
        END