New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
obs_sort.F90 in branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/OBS – NEMO

source: branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/OBS/obs_sort.F90 @ 5662

Last change on this file since 5662 was 5662, checked in by dancopsey, 9 years ago

Removed SVN keywords.

File size: 4.5 KB
Line 
1MODULE obs_sort
2   !!=====================================================================
3   !!                       ***  MODULE obs_sort  ***
4   !! Observation diagnostics: Various tools for sorting etc.
5   !!=====================================================================
6
7   !!----------------------------------------------------------------------
8   !!   sort_dp_indx : Get indicies for ascending order for a double prec. array
9   !!   index_sort   : Get indicies for ascending order for a double prec. array
10   !!---------------------------------------------------------------------
11   !! * Modules used
12   USE par_kind, ONLY : & ! Precision variables
13      & dp
14 
15   IMPLICIT NONE
16
17   !! * Routine accessibility
18   PRIVATE index_sort    ! Get indicies for ascending order for a double prec. array
19   
20   PUBLIC sort_dp_indx   ! Get indicies for ascending order for a double prec. array
21 
22   !!----------------------------------------------------------------------
23   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
24   !! $Id$
25   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
26   !!----------------------------------------------------------------------
27
28CONTAINS
29
30   SUBROUTINE sort_dp_indx( kvals, pvals, kindx )
31      !!----------------------------------------------------------------------
32      !!                    ***  ROUTINE sort_dp_indx  ***
33      !!         
34      !! ** Purpose : Get indicies for ascending order for a double precision array
35      !!
36      !! ** Method  : Call index_sort routine
37      !!
38      !! ** Action  :
39      !!
40      !! History :
41      !!        !  06-05  (K. Mogensen)  Original code
42      !!        !  06-10  (A. Weaver) Cleaning
43      !!----------------------------------------------------------------------
44
45      !! * Arguments
46      INTEGER, INTENT(IN) :: kvals     ! Number of elements to be sorted
47      REAL(KIND=dp), DIMENSION(kvals), INTENT(IN) :: &
48         & pvals            ! Array to be sorted
49      INTEGER, DIMENSION(kvals), INTENT(OUT) ::  &
50         & kindx            ! Indices for ordering of array
51
52      !! * Local declarations
53
54      !-----------------------------------------------------------------------
55      ! Call qsort routine
56      !-----------------------------------------------------------------------
57      IF (kvals>=1) THEN
58
59         CALL index_sort( pvals, kindx, kvals )
60
61      ENDIF
62
63   END SUBROUTINE sort_dp_indx
64
65   SUBROUTINE index_sort( pval, kindx, kvals )
66      !!----------------------------------------------------------------------
67      !!                    ***  ROUTINE index_sort  ***
68      !!         
69      !! ** Purpose : Get indicies for ascending order for a double precision array
70      !!
71      !! ** Method  : Heapsort
72      !!
73      !! ** Action  :
74      !!
75      !! References : http://en.wikipedia.org/wiki/Heapsort
76      !!
77      !! History :
78      !!        !  06-05  (K. Mogensen)  Original code
79      !!        !  06-10  (A. Weaver) Cleaning
80      !!----------------------------------------------------------------------
81
82      !! * Arguments
83      INTEGER, INTENT(IN) :: kvals         ! Number of values
84      REAL(KIND=dp), DIMENSION(kvals), INTENT(IN) :: &
85         & pval                            ! Array to be sorted
86      INTEGER, DIMENSION(kvals), INTENT(INOUT) :: &
87         & kindx                           ! Indicies for ordering
88
89      !! * Local declarations
90      INTEGER :: ji
91      INTEGER :: jj
92      INTEGER :: jt
93      INTEGER :: jn
94      INTEGER :: jparent
95      INTEGER :: jchild
96
97      DO ji = 1, kvals
98         kindx(ji) = ji
99      END DO
100     
101      ji = kvals/2 + 1
102      jn = kvals
103
104      main_loop : DO
105
106         IF ( ji > 1 ) THEN
107            ji = ji-1
108            jt = kindx(ji)
109         ELSE
110            jt = kindx(jn)
111            kindx(jn) = kindx(1)
112            jn = jn-1
113            IF ( jn <= 1 ) THEN
114               kindx(1) = jt
115               EXIT main_loop
116            ENDIF
117         ENDIF
118
119         jparent = ji
120         jchild =  2 * ji
121
122         inner_loop : DO
123
124            IF ( jchild > jn ) EXIT inner_loop
125            IF ( jchild < jn ) THEN
126               IF ( pval(kindx(jchild)) < pval(kindx(jchild+1)) ) THEN
127                 jchild = jchild+1
128               ENDIF
129            ENDIF
130            IF  ( pval(jt) < pval(kindx(jchild))) THEN
131               kindx(jparent) = kindx(jchild)
132               jparent = jchild
133               jchild  = jchild*2
134            ELSE
135               jchild = jn + 1 
136            ENDIF
137
138         END DO inner_loop
139
140         kindx(jparent) = jt
141
142      END DO main_loop
143     
144   END SUBROUTINE index_sort
145
146END MODULE obs_sort
147 
Note: See TracBrowser for help on using the repository browser.