这就是我最终要做的……我mrgrnk
从ORDERPACK中获取了代码,并根据我的目的进行了调整。下面的子程序findmap
似乎正在做我想做的事情。
module fndmap
use m_mrgrnk, only:mrgrnk
implicit none
contains
subroutine findmap(stkprm, stkmap )
! given 2-d real array stkprm, find a mapping described below:
!
! (identical records are assigned with same index)
! stkmap(i) == stkmap(j) iff stkprm(:,i) == stkprm(:,j)
! (order conserved)
! if i < j and stkmap(i) /= stkmap(j), then stkmap(i) < stkmap(j)
! (new index are contiguous)
! set(stkmap) == {1,2,..,maxval(stkmap)}
!
real,dimension(:,:),intent(in) :: stkprm
integer,dimension(:), intent(out) :: stkmap
integer, dimension(size(stkprm,2)) :: irngt
integer, dimension(size(stkprm,2)) :: iwork
integer :: nrec, i, j
nrec = size(stkprm,2)
! find rank of each record, duplicate records kept
call ar_mrgrnk(stkprm, irngt)
! construct iwork array, which has index of original array where the
! record are identical, and the index is youguest
i = 1
do while(i<=nrec)
do j=i+1,nrec
if (any(stkprm(:,irngt(i))/=stkprm(:,irngt(j)))) exit
enddo
iwork(irngt(i:j-1)) = minval(irngt(i:j-1))
i = j
enddo
! now construct the map, where stkmap(i) shows index of new array
! with duplicated record eliminated, original order kept
j = 0
do i=1,nrec
if (i==iwork(i)) then
j = j+1
stkmap(i) = j
else
stkmap(i) = stkmap(iwork(i))
endif
enddo
end subroutine
recursive subroutine ar_mrgrnk(xdont, irngt)
! behaves like mrgrnk of ORDERPACK, except that array is 2-d
! each row are ranked by first field, then second and so on
real, dimension(:,:), intent(in) :: xdont
integer, dimension(:), intent(out), target :: irngt
integer, dimension(size(xdont,2)) :: iwork
integer :: nfld,nrec
integer :: i, j
integer, dimension(:), pointer :: ipt
nfld=size(xdont,1)
nrec=size(xdont,2)
! rank by the first field
call mrgrnk(xdont(1,:), irngt)
! if there's only one field, it's done
if (nfld==1) return
! examine the rank to see if multiple record has identical
! values for the first field
i = 1
do while(i<=nrec)
do j=i+1,nrec
if (xdont(1,irngt(i))/=xdont(1,irngt(j))) exit
enddo
! if one-to-one, do nothing
if (j-1>i) then
! if many-to-one,
! gather those many, and rank them
call ar_mrgrnk(xdont(2:,irngt(i:j-1)),iwork)
! rearrange my rank based on those fields to the right
ipt => irngt(i:j-1)
ipt = ipt(iwork(1:j-i))
endif
i = j
enddo
if(associated(ipt)) nullify(ipt)
end subroutine
end module