module tracking
use precision
implicit none


public tag, partout
private

integer ntrackp(3)
real(lk),dimension(:,:),allocatable :: ptrackedi, ptrackede, ptrackedfki
save
contains
!=======================================================


subroutine tag
    use global_parameters, only: nparami, nparame, nparamfki
    use particle_array, only: mimax, memax, mfkimax
    use read_in_parameters, only: ifullk,irun
  
  ! Each MPI process has its own "ptracked" array to accumulate the tracked
  ! particles that may or may not reside in its subdomain.
  ! The vector "ntrackp" contains the number of tracked particles currently residing on a MPI process
  
    ntrackp=0
    if(ifullk==0)then
      allocate(ptrackedi(nparami,1)) !!array to store ion tracked particles
      ptrackedi=0.0_lk
    else
      allocate(ptrackedfki(nparamfki,1))  !!array to store tracked electrons
      ptrackedfki=0.0_lk
    endif
    allocate(ptrackede(nparame,1))  !!array to store tracked electrons
    ptrackede=0.0_lk
  
  ! If irun==0 and if particle tracking is "on", tag each particle with a unique number
  if(irun==0)  call tag_particles
  
end subroutine tag
  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  
subroutine partout
    call locate_tracked_particles
    call write_tracked_particles
end subroutine partout
  
!========================================
  !! originally written by s. either, modified by Y. Xiao April 2008
subroutine tag_particles
  
!========================================
  
    use global_parameters, only: mype,mi,me,nparami,nparame,mfki,nparamfki
    use particle_array, only: zion, zelectron, zion0, zelectron0,zfki,zfki0
    use read_in_parameters, only: ifullk

    integer :: m,np
  ! We keep track of the particles by tagging them with a unique number.
  ! We add an extra element to the particle array, which holds the
  ! particle tag, i.e. just a number.
  ! The input parameter "nptrack" is the total number of particles that
  ! we track.
  
  ! We tag each particle with a unique number that will be carried along
  ! with it as it moves from one processor/domain to another. To facilitate
  ! the search for the tracked particles, we give the non-tracked particles
  ! a negative value.
  ! CAVEAT: We are storing the tag in a floating point array, which means that
  !         for a large number of particles, the truncation error associated
  !         with the precision may make some particles have the same tag.
  !         This is particularly true in single precision where the maximum
  !         number of particles that can be distinguished with tags differing
  !         by unity is 2**24-1 = 16,777,215 (16.7 million). The most important
  !         is to distinguish the tracked particles between them. The other
  !         particles just need to have a negative tag. In order to minimize
  !         the effect of the precision truncation, we retag the tracked
  !         particles with positive numbers starting at 1.
  
    do m=1,me
        zelectron(nparame-1,m)=-real(m,lk) !!!-real(m+mype*mi)
        zelectron(nparame,m)=real(mype+1,lk)
        zelectron0(nparame-1,m)=-real(m,lk)
        zelectron0(nparame,m)=real(mype+1,lk)
    enddo
    if(ifullk==0)then
      do m=1,mi
        zion(nparami-1,m)=-real(m,lk) !!!-real(m+mype*mi)
        zion(nparami,m)=real(mype+1,lk)
        zion0(nparami-1,m)=-real(m,lk)
        zion0(nparami,m)=real(mype+1,lk)
      enddo
    else
      do m=1,mfki
        zfki(nparamfki-1,m)=-real(m,lk)
        zfki(nparamfki,m)=real(mype+1,lk)
        zfki0(nparamfki-1,m)=-real(m,lk)
        zfki0(nparamfki,m)=real(mype+1,lk)
      enddo
    endif

    np=0
    if(mype==0)then
      m=1
        np=np+1
        zelectron(nparame-1,m)=real(np,lk)
        zelectron0(nparame-1,m)=real(np,lk)
    endif
  
    if(mype==0)write(*,*)'np=',np
  !!  if(mype==0)write(*,*)'zelectron(1,1:20)=',zelectron(1,1:20),'zelectron(nparame,1)=',zelectron(nparame,1)
  !! we only condiser those particles originally satify a certain condition
  ! namely within an annulous ain<r<aout
  ! We divide the number of tracked particles equally between all processors
  ! as much as possible. If nptrack is not a multiple of numberpe, we add
  ! 1 particle to some of the the processors until we get to nptrack. We also
  ! start at mi/2 so that the lower numbers fall around r=0.5 (see subroutine 
  ! load for the initial particle distribution in r).
  
    np=0
   ! this part is for tagging a single particle 
   ! if(mype==0)then
   ! np=np+1
   ! zion(nparami-1,1)=real(np,wp)
          ! zion0(nparami-1,1)=real(np,wp)
   ! endif
  
  ! use this part for selecting a group of particles
    if(ifullk==0)then
      if(mype==0)then
        m=1
        np=np+1
        zion(nparami-1,m)=real(np,lk)
        zion0(nparami-1,m)=real(np,lk)
      endif
    else
      if(mype==0)then
        m=mfki/2
        np=np+1
        zfki(nparamfki-1,m)=real(np,lk)
        zfki0(nparamfki-1,m)=real(np,lk)
      endif
    endif
  
  !  if(mype==0)write(78,*)'np=',np,'  npp=',npp,'  zion(nparami-1,mi/2)=',zion(nparami-1,mi/2),&
  !                        '  zion(nparami-1,1)=',zion(nparami-1,1)
  
  !  if(mype==0)then
  !  ! On the master process (mype=0), we pick "nptrack" particles that will
  !  ! be followed at every time step. To facilitate the search of those
  !  ! particles among all the others in the zion arrays of each processor,
  !  ! we give them a positive number from 1 to nptrack. The particles are
  !  ! picked around r=0.5, which is mi/2.
  !    do m=(mi-nptrack)/2,(mi+nptrack)/2-1
  !       zion(nparami-1,m)=-zion(nparami-1,m)
  !       zion0(nparami-1,m)=-zion(nparami-1,m)
  !    enddo
  !  endif

  
end subroutine tag_particles
  
  !========================================
  
subroutine locate_tracked_particles
  
  !========================================
  
    use global_parameters, only: mi, nparami, me, nparame, mfki,nparamfki
    use particle_array, only: zion, zelectron, zfki
    use read_in_parameters, only: ifullk, kinetice
    integer :: m,npp
  
  ! Check if tracked particles are located on this processor. Particles that
  ! are tracked at every time step have a positive zion(nparami-1,m) value. All the
  ! others have a negative value.
    ntrackp=0
    if(ifullk==0)then
    ptrackedi=0.0
    else
    ptrackedfki=0.0
    endif
    npp=0
    if(ifullk==0)then 
        do m=1,mi
            if(zion(nparami-1,m)>0.0)then
                npp=npp+1
                ptrackedi(1:nparami,npp)=zion(1:nparami,m)
            endif
        enddo
    else
        do m=1,mfki
            if(zfki(nparamfki-1,m)>0.0)then
                npp=npp+1
                ptrackedfki(1:nparamfki,npp)=zfki(1:nparamfki,m)
            endif
        enddo
    endif
    ntrackp(1)=npp
    npp=0
  
    ptrackede=0.0
    do m=1,me
        if(zelectron(nparame-1,m)>0.0)then
            npp=npp+1
            ptrackede(1:nparame,npp)=zelectron(1:nparame,m)
        endif
    enddo
    ntrackp(2)=npp
    if(kinetice==0) ntrackp(2)=0
  !!  if(mype==0)write(*,*)'me=',me,'np=',npp,'zelectron(1,1:6)=',zelectron(nparame-1,1:6),&
  !!                        'zelectron(nparame,1)= ',zelectron(nparame,1)
  
end subroutine locate_tracked_particles
  
  !========================================
  
subroutine write_tracked_particles
  
  !========================================
  
    use global_parameters, only: istep, nspecies, mype, nparami, nparame, nparamfki
    use read_in_parameters, only: irun, mstep, ndiag, ifullk
    integer :: j
    character(len=30) :: cdum
    write(cdum,'("trackp_dir/TRACKP.",i5.5)')mype
    cdum=trim(cdum)//char(0)
    if(irun==0 .and. istep==ndiag)then
        open(57,file=cdum,status='replace',position='append')
    else
        open(57,file=cdum,status='old',position='append')
    endif
    write(57,*)istep+irun*mstep
    write(57,*)ntrackp(1:nspecies)
    if(ifullk==0)then
        do j=1,ntrackp(1)
            write(57,'(12(E17.9)/)')ptrackedi(1:nparami,j)
        enddo
    else
        do j=1,ntrackp(1)
            write(57,'(15(E17.9)/)')ptrackedfki(1:nparamfki,j)
        enddo
    endif
    do j=1,ntrackp(2)
        write(57,'(12(E17.9)/)')ptrackede(1:nparame,j)
    enddo

    close(57)
end subroutine write_tracked_particles
  
end module tracking
