*IF DEF,A13_1B                                                             QTPOS1B.2      
C ******************************COPYRIGHT******************************    QTPOS1B.3      
C (c) CROWN COPYRIGHT 1996, METEOROLOGICAL OFFICE, All Rights Reserved.    QTPOS1B.4      
C                                                                          QTPOS1B.5      
C Use, duplication or disclosure of this code is subject to the            QTPOS1B.6      
C restrictions as set forth in the contract.                               QTPOS1B.7      
C                                                                          QTPOS1B.8      
C                Meteorological Office                                     QTPOS1B.9      
C                London Road                                               QTPOS1B.10     
C                BRACKNELL                                                 QTPOS1B.11     
C                Berkshire UK                                              QTPOS1B.12     
C                RG12 2SZ                                                  QTPOS1B.13     
C                                                                          QTPOS1B.14     
C If no contract has been raised with this copy of the code, the use,      QTPOS1B.15     
C duplication or disclosure of it is strictly prohibited.  Permission      QTPOS1B.16     
C to do so must first be obtained in writing from the Head of Numerical    QTPOS1B.17     
C Modelling at the above address.                                          QTPOS1B.18     
C ******************************COPYRIGHT******************************    QTPOS1B.19     
C                                                                          QTPOS1B.20     
!+ Interface routine for QTPOS allowing both normal and MPP code to use    QTPOS1B.21     
!+ current algorithms                                                      QTPOS1B.22     
!                                                                          QTPOS1B.23     
! Subroutine Interface                                                     QTPOS1B.24     

      SUBROUTINE QT_POS_CTL                                                 1,12QTPOS1B.25     
     &                     (QT,RS_SQUARED_DELTAP,ROW_LENGTH,               QTPOS1B.26     
     &                      P_FIELD,Q_LEVELS,                              QTPOS1B.27     
*CALL ARGFLDPT                                                             QTPOS1B.28     
     &                      ERROR_CODE,ERROR_MESSAGE,                      QTPOS1B.29     
     &                      COS_P_LATITUDE,SEC_P_LATITUDE,                 QTPOS1B.30     
     &                      L_NEG_QT,L_QT_POS_LOCAL,DT,SF_QTFIX,QTFIX)     QTPOS1B.31     
                                                                           QTPOS1B.32     
      IMPLICIT NONE                                                        QTPOS1B.33     
                                                                           QTPOS1B.34     
!                                                                          QTPOS1B.35     
! Description:                                                             QTPOS1B.36     
! QT fields are mass weighted, and then                                    QTPOS1B.37     
! if MPP *DEF is set, the field is gathered a level to a processor, and    QTPOS1B.38     
! QT_POS is called by each processor with its levels. Otherwise (ie. non   QTPOS1B.39     
! MPP), this routine just passed the mass-weighted QT field straight       QTPOS1B.40     
! through to QTPOS. Finally, the mass weighting is removed, and            QTPOS1B.41     
! the QT_FIX diagnostic is calculated (if required)                        QTPOS1B.42     
!                                                                          QTPOS1B.43     
! Current code owner (QT_POS_CTL only) : Paul Burton                       QTPOS1B.44     
!                                                                          QTPOS1B.45     
! History                                                                  QTPOS1B.46     
!  Model    Date      Modification history:                                QTPOS1B.47     
!   4.2   28/10/96  New deck for HADCM2-specific section A13_1B,           QTPOS1B.48     
!                   as QTPOS1A but with inconsistent 'old' type            QTPOS1B.49     
!                   of polar weights.  T.Johns                             QTPOS1B.50     
!    4.3    24/04/97  Improved error trapping for MPP code                 GPB5F403.73     
!                     Set QTFIX to safe values         P.Burton            GPB5F403.74     
!    4.4    22/09/97  Correct missetting of local variables  T Johns       ATJ0F404.1      
!    4.5    29/09/98  T3E only: remove negative zeros from QT              GPB0F405.196    
!                                                    P.Burton              GPB0F405.197    
!                                                                          QTPOS1B.51     
! Subroutine Arguments:                                                    QTPOS1B.52     
                                                                           QTPOS1B.53     
      INTEGER                                                              QTPOS1B.54     
     &  P_FIELD                  ! IN : Size of fields on P grid           QTPOS1B.55     
     &, ROW_LENGTH               ! IN : number of points per row           QTPOS1B.56     
     &, Q_LEVELS                 ! IN : number of moist levels             QTPOS1B.57     
!                                                                          QTPOS1B.58     
     &, ERROR_CODE               ! INOUT : 0 on entry                      QTPOS1B.59     
!                                !         1 on exit if error detected     QTPOS1B.60     
                                                                           QTPOS1B.61     
! All TYPFLDPT arguments are intent IN                                     QTPOS1B.62     
*CALL TYPFLDPT                                                             QTPOS1B.63     
                                                                           QTPOS1B.64     
      CHARACTER*(80)                                                       QTPOS1B.65     
     &  ERROR_MESSAGE            ! OUT : Error message                     QTPOS1B.66     
                                                                           QTPOS1B.67     
      REAL                                                                 QTPOS1B.68     
     &  QT(P_FIELD,Q_LEVELS)     ! INOUT : QT field                        QTPOS1B.69     
     &, QTFIX(P_FIELD,Q_LEVELS) ! OUT : QT fix diagnostic from             QTPOS1B.70     
!                                !       resetting and redistributing      QTPOS1B.71     
!                                !       negative QT values                QTPOS1B.72     
     &, RS_SQUARED_DELTAP(P_FIELD,Q_LEVELS)                                QTPOS1B.73     
!                                ! IN : RS*RS*DELTA_P                      QTPOS1B.74     
     &, COS_P_LATITUDE(P_FIELD)  ! IN : COS(LATITUDE) at P points          QTPOS1B.75     
     &, SEC_P_LATITUDE(P_FIELD)  ! IN : SEC(LATITUDE) at P points          QTPOS1B.76     
     &, DT                       ! IN : dynamics timestep                  QTPOS1B.77     
                                                                           QTPOS1B.78     
      LOGICAL                                                              QTPOS1B.79     
     &  L_NEG_QT                 ! IN : Switch to continue in event of     QTPOS1B.80     
!                                !      excessive negative QT              QTPOS1B.81     
!                                !      (non-conservative)                 QTPOS1B.82     
     &, L_QT_POS_LOCAL           ! IN : perform local correction           QTPOS1B.83     
!                                !      (method 2)                         QTPOS1B.84     
     &, SF_QTFIX                 ! IN : STASHflag for output of QT fix     QTPOS1B.85     
!                                !      diagnostic                         QTPOS1B.86     
                                                                           QTPOS1B.87     
! Local variables                                                          QTPOS1B.88     
                                                                           QTPOS1B.89     
      INTEGER                                                              QTPOS1B.90     
     &  K,I                      ! loop variables                          QTPOS1B.91     
     &, I_start                  ! start address for QTPOS calcs           QTPOS1B.92     
     &, I_end                    ! end address of QTPOS calcs              QTPOS1B.93     
*IF DEF,MPP                                                                QTPOS1B.94     
     &, I_start_local            ! local start address for QTPOS calcs     QTPOS1B.95     
     &, I_end_local              ! local end address for QTPOS calcs       QTPOS1B.96     
     &, info                     ! return code from comms.                 QTPOS1B.97     
     &, MAP(Q_LEVELS)            ! processor number for level              QTPOS1B.98     
     &, N_LEVS_ON_PROC(0:N_PROCS-1) ! number of levels on each processor   QTPOS1B.99     
     &, int_FAILURE(Q_LEVELS)    ! version of LOGICAL FAILURE array,       QTPOS1B.100    
!                                ! but using integers                      QTPOS1B.101    
*ENDIF                                                                     QTPOS1B.102    
                                                                           QTPOS1B.103    
      REAL                                                                 QTPOS1B.104    
     &  FACTOR(Q_LEVELS)         ! scaling factor for revised QT           QTPOS1B.105    
*IF DEF,GLOBAL                                                             QTPOS1B.106    
     &, LATITUDE_STEP_INVERSE    ! Computed 1/(DELTA_LAMBDA)               QTPOS1B.107    
     &, COS_P_LAT_NP             ! COS(latitude) at North Pole             QTPOS1B.108    
     &, COS_P_LAT_SP             ! COS(latitude) at South Pole             QTPOS1B.109    
     &, SEC_P_LAT_NP             ! SEC(latitude) at North Pole             QTPOS1B.110    
     &, SEC_P_LAT_SP             ! SEC(latitude) at South Pole             QTPOS1B.111    
! See note below concerning rotated global grids                           QTPOS1B.112    
*ENDIF                                                                     QTPOS1B.113    
*IF DEF,MPP                                                                QTPOS1B.114    
     &, local_FACTOR((Q_LEVELS/N_PROCS)+1)  ! scaling factors for the      QTPOS1B.115    
!                                          ! levels I process              QTPOS1B.116    
     &,  global_QT_data(GLOBAL_P_FIELD,(Q_LEVELS/N_PROCS)+1)               QTPOS1B.117    
! Array containing the QT levels to be processed on this processor         QTPOS1B.118    
     &,  NP_VALUE(Q_LEVELS,2)    ! value of QT + QTFIX at North Pole       QTPOS1B.119    
     &,  SP_VALUE(Q_LEVELS,2)    ! value of QT + QTFIX at South Pole       QTPOS1B.120    
*ENDIF                                                                     QTPOS1B.121    
                                                                           QTPOS1B.122    
      LOGICAL                                                              QTPOS1B.123    
     &  FAILURE(Q_LEVELS)  ! Indicates if a particular level has           QTPOS1B.124    
!                          ! failed to be QT_POSd.                         QTPOS1B.125    
*IF DEF,MPP                                                                QTPOS1B.126    
     &, local_FAILURE((Q_LEVELS/N_PROCS)+1)  ! FAILURE array for           QTPOS1B.127    
!                                           ! levels I process             QTPOS1B.128    
*ENDIF                                                                     QTPOS1B.129    
                                                                           QTPOS1B.130    
*CALL C_A                                                                  QTPOS1B.131    
*CALL C_G                                                                  QTPOS1B.132    
*CALL C_PI                                                                 QTPOS1B.133    
                                                                           QTPOS1B.134    
!---------------------------------------------------------------------     QTPOS1B.135    
! NOTE : The MPP code assumes, that for global models, the value of        QTPOS1B.136    
!        latitude remains constant over a row (in particular the polar     QTPOS1B.137    
!        rows). This implies the code will not give the correct answers    QTPOS1B.138    
!        if a rotated global grid is used.                                 QTPOS1B.139    
                                                                           QTPOS1B.140    
*IF DEF,GLOBAL                                                             QTPOS1B.141    
! Recalculate LATITUDE_STEP_INVERSE to avoid changing argument list        QTPOS1B.142    
! (HADCM2-specific code, needs *CALL C_PI)                                 QTPOS1B.143    
                                                                           QTPOS1B.144    
      LATITUDE_STEP_INVERSE=RECIP_PI_OVER_180 *                            ATJ0F404.2      
     &           REAL((GLOBAL_P_FIELD/GLOBAL_ROW_LENGTH)-1) / REAL(180)    ATJ0F404.3      
                                                                           QTPOS1B.147    
      IF (L_QT_POS_LOCAL) THEN                                             QTPOS1B.148    
        I_start=1                                                          QTPOS1B.149    
        I_end=GLOBAL_P_FIELD                                               QTPOS1B.150    
*IF DEF,MPP                                                                QTPOS1B.151    
        I_start_local=FIRST_FLD_PT                                         QTPOS1B.152    
        I_end_local=LAST_P_FLD_PT                                          QTPOS1B.153    
*ENDIF                                                                     QTPOS1B.154    
      ELSE ! global correction scheme                                      QTPOS1B.155    
        I_start=GLOBAL_ROW_LENGTH                                          QTPOS1B.156    
        I_end=GLOBAL_P_FIELD+1-GLOBAL_ROW_LENGTH                           QTPOS1B.157    
*IF DEF,MPP                                                                QTPOS1B.158    
        IF (at_top_of_LPG) THEN                                            QTPOS1B.159    
          I_start_local=TOP_ROW_START+LAST_ROW_PT-1                        QTPOS1B.160    
        ELSE                                                               QTPOS1B.161    
          I_start_local=FIRST_FLD_PT                                       QTPOS1B.162    
        ENDIF                                                              QTPOS1B.163    
        IF (at_base_of_LPG) THEN                                           QTPOS1B.164    
          I_end_local=P_BOT_ROW_START+FIRST_ROW_PT-1                       QTPOS1B.165    
        ELSE                                                               QTPOS1B.166    
          I_end_local=LAST_P_FLD_PT                                        QTPOS1B.167    
        ENDIF                                                              QTPOS1B.168    
*ENDIF                                                                     QTPOS1B.169    
                                                                           QTPOS1B.170    
*IF DEF,MPP                                                                QTPOS1B.171    
        IF (at_top_of_LPG .AND. at_right_of_LPG) THEN                      QTPOS1B.172    
*ENDIF                                                                     QTPOS1B.173    
          COS_P_LAT_NP=COS_P_LATITUDE(1)                                   QTPOS1B.174    
          SEC_P_LAT_NP=SEC_P_LATITUDE(1)                                   QTPOS1B.175    
          COS_P_LATITUDE(TOP_ROW_START+LAST_ROW_PT-1)=                     QTPOS1B.176    
     &      .125*GLOBAL_ROW_LENGTH/LATITUDE_STEP_INVERSE                   QTPOS1B.177    
          SEC_P_LATITUDE(TOP_ROW_START+LAST_ROW_PT-1)=                     QTPOS1B.178    
     &      1.0/COS_P_LATITUDE(TOP_ROW_START+LAST_ROW_PT-1)                QTPOS1B.179    
*IF DEF,MPP                                                                QTPOS1B.180    
        ENDIF                                                              QTPOS1B.181    
                                                                           QTPOS1B.182    
        IF (at_base_of_LPG .AND. at_left_of_LPG) THEN                      QTPOS1B.183    
*ENDIF                                                                     QTPOS1B.184    
       COS_P_LAT_SP=COS_P_LATITUDE(P_BOT_ROW_START+FIRST_ROW_PT-1)         ATJ0F404.4      
       SEC_P_LAT_SP=SEC_P_LATITUDE(P_BOT_ROW_START+FIRST_ROW_PT-1)         ATJ0F404.5      
          COS_P_LATITUDE(P_BOT_ROW_START+FIRST_ROW_PT-1)=                  QTPOS1B.187    
     &      .125*GLOBAL_ROW_LENGTH/LATITUDE_STEP_INVERSE                   QTPOS1B.188    
          SEC_P_LATITUDE(P_BOT_ROW_START+FIRST_ROW_PT-1)=                  QTPOS1B.189    
     &      1.0/COS_P_LATITUDE(P_BOT_ROW_START+FIRST_ROW_PT-1)             QTPOS1B.190    
*IF DEF,MPP                                                                QTPOS1B.191    
        ENDIF                                                              QTPOS1B.192    
*ENDIF                                                                     QTPOS1B.193    
                                                                           QTPOS1B.194    
      ENDIF                                                                QTPOS1B.195    
                                                                           QTPOS1B.196    
*ELSE                                                                      QTPOS1B.197    
      I_start=1                                                            QTPOS1B.198    
      I_end=GLOBAL_P_FIELD                                                 QTPOS1B.199    
*IF DEF,MPP                                                                QTPOS1B.200    
      I_start_local=FIRST_FLD_PT                                           QTPOS1B.201    
      I_end_local=LAST_P_FLD_PT                                            QTPOS1B.202    
*ENDIF                                                                     QTPOS1B.203    
*ENDIF                                                                     QTPOS1B.204    
                                                                           QTPOS1B.205    
      DO K=1,Q_LEVELS                                                      QTPOS1B.206    
*IF -DEF,MPP                                                               QTPOS1B.207    
        DO I=I_start,I_end                                                 QTPOS1B.208    
*ELSE                                                                      QTPOS1B.209    
        DO I=I_start_local,I_end_local                                     QTPOS1B.210    
*ENDIF                                                                     QTPOS1B.211    
          QT(I,K)=QT(I,K)*RS_SQUARED_DELTAP(I,K)*COS_P_LATITUDE(I)         QTPOS1B.212    
        ENDDO                                                              QTPOS1B.213    
                                                                           QTPOS1B.214    
        IF (SF_QTFIX) THEN                                                 QTPOS1B.215    
          DO I=FIRST_VALID_PT,LAST_P_VALID_PT                              QTPOS1B.216    
            QTFIX(I,K)=QT(I,K)                                             QTPOS1B.217    
          ENDDO                                                            QTPOS1B.218    
          DO I=1,FIRST_VALID_PT-1                                          GPB5F403.75     
            QTFIX(I,K)=QT(FIRST_VALID_PT,K)                                GPB5F403.76     
          ENDDO                                                            GPB5F403.77     
          DO I=LAST_P_VALID_PT+1,P_FIELD                                   GPB5F403.78     
            QTFIX(I,K)=QT(LAST_P_VALID_PT,K)                               GPB5F403.79     
          ENDDO                                                            GPB5F403.80     
        ENDIF                                                              QTPOS1B.219    
                                                                           QTPOS1B.220    
*IF DEF,MPP                                                                QTPOS1B.221    
! Caclulate the mapping of which processor each level will go to           QTPOS1B.222    
        MAP(K)=MOD((K-1),N_PROCS)  ! assumes first PE is PE 0              QTPOS1B.223    
                                                                           QTPOS1B.224    
*ENDIF                                                                     QTPOS1B.225    
      ENDDO                                                                QTPOS1B.226    
                                                                           QTPOS1B.227    
*IF DEF,MPP                                                                QTPOS1B.228    
! Distribute QT over the processors                                        QTPOS1B.229    
      DO K=0,N_PROCS-1                                                     QTPOS1B.230    
        N_LEVS_ON_PROC(K)=0                                                QTPOS1B.231    
      ENDDO                                                                QTPOS1B.232    
                                                                           QTPOS1B.233    
      DO K=1,Q_LEVELS                                                      QTPOS1B.234    
        N_LEVS_ON_PROC(MAP(K))=N_LEVS_ON_PROC(MAP(K))+1                    QTPOS1B.235    
                                                                           QTPOS1B.236    
        CALL GATHER_FIELD(QT(1,K),                                         QTPOS1B.237    
     &                    global_QT_data(1,N_LEVS_ON_PROC(MAP(K))),        QTPOS1B.238    
     &                    ROW_LENGTH,tot_P_ROWS,                           QTPOS1B.239    
     &                    GLOBAL_ROW_LENGTH,                               QTPOS1B.240    
     &                    GLOBAL_P_FIELD/GLOBAL_ROW_LENGTH,                QTPOS1B.241    
     &                    MAP(K),GC_ALL_GROUP,info)                        QTPOS1B.242    
                                                                           QTPOS1B.243    
      ENDDO                                                                QTPOS1B.244    
                                                                           QTPOS1B.245    
! and call QT_POS with the whole levels on this processor                  QTPOS1B.246    
                                                                           QTPOS1B.247    
      DO K=1,N_LEVS_ON_PROC(MY_PROC_ID)                                    QTPOS1B.248    
        local_FAILURE(K)=.FALSE.                                           QTPOS1B.249    
      ENDDO                                                                QTPOS1B.250    
                                                                           QTPOS1B.251    
      CALL QT_POS(global_QT_data,global_ROW_LENGTH,global_P_FIELD,         QTPOS1B.252    
     &            I_start,I_end,                                           QTPOS1B.253    
     &            N_LEVS_ON_PROC(MY_PROC_ID),                              QTPOS1B.254    
     &            ERROR_CODE,ERROR_MESSAGE,                                QTPOS1B.255    
     &            local_FACTOR,local_FAILURE,                              QTPOS1B.256    
     &            L_NEG_QT,L_QT_POS_LOCAL)                                 QTPOS1B.257    
                                                                           QTPOS1B.258    
! Set ERROR_CODE and ERROR_MESSAGE if any pe has set ERROR_CODE            GPB5F403.81     
                                                                           GPB5F403.82     
      CALL GC_IMAX(1,n_procs,info,ERROR_CODE)                              GPB5F403.83     
                                                                           GPB5F403.84     
      IF (ERROR_CODE .NE. 0) THEN                                          GPB5F403.85     
        ERROR_MESSAGE=                                                     GPB5F403.86     
     & 'QT_POS   : MASS-WEIGHTED QT SUMMED OVER A LEVEL WAS NEGATIVE.'     GPB5F403.87     
      ENDIF                                                                GPB5F403.88     
! and gather back levels to original processors                            QTPOS1B.259    
                                                                           QTPOS1B.260    
      DO K=0,N_PROCS-1                                                     QTPOS1B.261    
        N_LEVS_ON_PROC(K)=0                                                QTPOS1B.262    
      ENDDO                                                                QTPOS1B.263    
                                                                           QTPOS1B.264    
      DO K=1,Q_LEVELS                                                      QTPOS1B.265    
        N_LEVS_ON_PROC(MAP(K))=N_LEVS_ON_PROC(MAP(K))+1                    QTPOS1B.266    
                                                                           QTPOS1B.267    
! collect the FAILURE and FACTOR array elements for this level             QTPOS1B.268    
                                                                           QTPOS1B.269    
        IF (MY_PROC_ID .EQ. MAP(K)) THEN                                   QTPOS1B.270    
!         I processed this level                                           QTPOS1B.271    
          IF (local_FAILURE(N_LEVS_ON_PROC(MY_PROC_ID))) THEN              QTPOS1B.272    
            int_FAILURE(K)=1                                               QTPOS1B.273    
          ELSE                                                             QTPOS1B.274    
            int_FAILURE(K)=0                                               QTPOS1B.275    
          ENDIF                                                            QTPOS1B.276    
          FACTOR(K)=local_FACTOR(N_LEVS_ON_PROC(MY_PROC_ID))               QTPOS1B.277    
        ELSE  ! I didn't process this level                                QTPOS1B.278    
          int_FAILURE(K)=0                                                 QTPOS1B.279    
          FACTOR(K)=0.0                                                    QTPOS1B.280    
        ENDIF                                                              QTPOS1B.281    
                                                                           QTPOS1B.282    
! and scatter this level back to it's original home                        QTPOS1B.283    
                                                                           QTPOS1B.284    
        CALL SCATTER_FIELD(QT(1,K),                                        QTPOS1B.285    
     &                    global_QT_data(1,N_LEVS_ON_PROC(MAP(K))),        QTPOS1B.286    
     &                    ROW_LENGTH,tot_P_ROWS,                           QTPOS1B.287    
     &                    GLOBAL_ROW_LENGTH,                               QTPOS1B.288    
     &                    GLOBAL_P_FIELD/GLOBAL_ROW_LENGTH,                QTPOS1B.289    
     &                    MAP(K),GC_ALL_GROUP,info)                        QTPOS1B.290    
                                                                           QTPOS1B.291    
      ENDDO                                                                QTPOS1B.292    
                                                                           QTPOS1B.293    
! and distribute the FAILURE and FACTOR arrays                             QTPOS1B.294    
                                                                           QTPOS1B.295    
      CALL GC_ISUM(Q_LEVELS,N_PROCS,info,int_FAILURE)                      QTPOS1B.296    
      CALL GC_RSUM(Q_LEVELS,N_PROCS,info,FACTOR)                           QTPOS1B.297    
                                                                           QTPOS1B.298    
      DO K=1,Q_LEVELS                                                      QTPOS1B.299    
        IF (int_FAILURE(K) .EQ. 0) THEN                                    QTPOS1B.300    
          FAILURE(K)=.FALSE.                                               QTPOS1B.301    
        ELSE                                                               QTPOS1B.302    
          FAILURE(K)=.TRUE.                                                QTPOS1B.303    
        ENDIF                                                              QTPOS1B.304    
      ENDDO                                                                QTPOS1B.305    
                                                                           QTPOS1B.306    
*ELSE                                                                      QTPOS1B.307    
! Call QT_POS                                                              QTPOS1B.308    
                                                                           QTPOS1B.309    
      CALL QT_POS(QT,ROW_LENGTH,P_FIELD,I_start,I_end,Q_LEVELS,            QTPOS1B.310    
     &            ERROR_CODE,ERROR_MESSAGE,FACTOR,FAILURE,                 QTPOS1B.311    
     &            L_NEG_QT,L_QT_POS_LOCAL)                                 QTPOS1B.312    
                                                                           QTPOS1B.313    
*ENDIF                                                                     QTPOS1B.314    
                                                                           QTPOS1B.315    
! Now remove the mass weighting, and calculate QT_FIX diagnostic           QTPOS1B.316    
                                                                           QTPOS1B.317    
      DO K=1,Q_LEVELS                                                      QTPOS1B.318    
        IF (FAILURE(K)) THEN                                               QTPOS1B.319    
          IF (SF_QTFIX) THEN                                               QTPOS1B.320    
*IF -DEF,MPP                                                               QTPOS1B.321    
            DO I=I_start,I_end                                             QTPOS1B.322    
*ELSE                                                                      QTPOS1B.323    
            DO I=I_start_local,I_end_local                                 QTPOS1B.324    
*ENDIF                                                                     QTPOS1B.325    
              QTFIX(I,K) = (QT(I,K) * FACTOR(K) - QTFIX(I,K))*             QTPOS1B.326    
     &                     SEC_P_LATITUDE(I)/(A*A*G*DT)                    QTPOS1B.327    
            ENDDO                                                          QTPOS1B.328    
          ENDIF                                                            QTPOS1B.329    
*IF -DEF,MPP                                                               QTPOS1B.330    
          DO I=I_start,I_end                                               QTPOS1B.331    
*ELSE                                                                      QTPOS1B.332    
          DO I=I_start_local,I_end_local                                   QTPOS1B.333    
*ENDIF                                                                     QTPOS1B.334    
            QT(I,K)=QT(I,K)*FACTOR(K)/RS_SQUARED_DELTAP(I,K)*              QTPOS1B.335    
     &              SEC_P_LATITUDE(I)                                      QTPOS1B.336    
          ENDDO                                                            QTPOS1B.337    
        ELSE  ! FAILURE(K) .EQ. .FALSE.                                    QTPOS1B.338    
          IF (SF_QTFIX) THEN                                               QTPOS1B.339    
*IF -DEF,MPP                                                               QTPOS1B.340    
            DO I=I_start,I_end                                             QTPOS1B.341    
*ELSE                                                                      QTPOS1B.342    
            DO I=I_start_local,I_end_local                                 QTPOS1B.343    
*ENDIF                                                                     QTPOS1B.344    
              QTFIX(I,K) = (QT(I,K) - QTFIX(I,K)) *                        QTPOS1B.345    
     &                     SEC_P_LATITUDE(I)/(A*A*G*DT)                    QTPOS1B.346    
            ENDDO                                                          QTPOS1B.347    
          ENDIF                                                            QTPOS1B.348    
*IF -DEF,MPP                                                               QTPOS1B.349    
          DO I=I_start,I_end                                               QTPOS1B.350    
*ELSE                                                                      QTPOS1B.351    
          DO I=I_start_local,I_end_local                                   QTPOS1B.352    
*ENDIF                                                                     QTPOS1B.353    
            QT(I,K)=QT(I,K)/RS_SQUARED_DELTAP(I,K)*                        QTPOS1B.354    
     &              SEC_P_LATITUDE(I)                                      QTPOS1B.355    
          ENDDO                                                            QTPOS1B.356    
        ENDIF                                                              QTPOS1B.357    
                                                                           QTPOS1B.358    
      ENDDO ! K : loop over levels                                         QTPOS1B.359    
                                                                           QTPOS1B.360    
*IF DEF,GLOBAL                                                             QTPOS1B.361    
      IF (.NOT. L_QT_POS_LOCAL) THEN                                       QTPOS1B.362    
                                                                           QTPOS1B.363    
! Set all points at poles equal to last(NP)/first(SP) point                QTPOS1B.364    
*IF -DEF,MPP                                                               QTPOS1B.365    
          DO K=1,Q_LEVELS                                                  QTPOS1B.366    
            DO I=1,ROW_LENGTH-1                                            QTPOS1B.367    
              QT(I,K) = QT(ROW_LENGTH,K)                                   QTPOS1B.368    
              QT(P_FIELD+1-I,K) = QT(P_FIELD+1-ROW_LENGTH,K)               QTPOS1B.369    
            ENDDO                                                          QTPOS1B.370    
                                                                           QTPOS1B.371    
            IF (SF_QTFIX) THEN                                             QTPOS1B.372    
              DO I=1,ROW_LENGTH-1                                          QTPOS1B.373    
                QTFIX(I,K) = QTFIX(ROW_LENGTH,K)                           QTPOS1B.374    
                QTFIX(P_FIELD+1-I,K) = QTFIX(P_FIELD+1-ROW_LENGTH,K)       QTPOS1B.375    
              ENDDO                                                        QTPOS1B.376    
            ENDIF                                                          QTPOS1B.377    
          ENDDO                                                            QTPOS1B.378    
*ELSE                                                                      QTPOS1B.379    
! Set up arrays containing NP values                                       QTPOS1B.380    
          IF (at_top_of_LPG .AND. at_right_of_LPG) THEN                    QTPOS1B.381    
            DO K=1,Q_LEVELS                                                QTPOS1B.382    
              NP_VALUE(K,1)=QT(I_start_local,K)                            QTPOS1B.383    
              NP_VALUE(K,2)=QTFIX(I_start_local,K)                         QTPOS1B.384    
            ENDDO                                                          QTPOS1B.385    
          ELSE                                                             QTPOS1B.386    
            DO K=1,Q_LEVELS                                                QTPOS1B.387    
              NP_VALUE(K,1)=0.0                                            QTPOS1B.388    
              NP_VALUE(K,2)=0.0                                            QTPOS1B.389    
            ENDDO                                                          QTPOS1B.390    
          ENDIF                                                            QTPOS1B.391    
                                                                           QTPOS1B.392    
! Distribute the NP value over processors on polar row                     QTPOS1B.393    
          IF (at_top_of_LPG) THEN                                          QTPOS1B.394    
            CALL GCG_RSUM(Q_LEVELS,GC_ROW_GROUP,info,NP_VALUE(1,1))        QTPOS1B.395    
            DO K=1,Q_LEVELS                                                QTPOS1B.396    
              DO I=TOP_ROW_START,TOP_ROW_START+ROW_LENGTH-1                QTPOS1B.397    
                QT(I,K)=NP_VALUE(K,1)                                      QTPOS1B.398    
              ENDDO                                                        QTPOS1B.399    
            ENDDO                                                          QTPOS1B.400    
            IF (SF_QTFIX) THEN                                             QTPOS1B.401    
              CALL GCG_RSUM(Q_LEVELS,GC_ROW_GROUP,info,NP_VALUE(1,2))      QTPOS1B.402    
              DO K=1,Q_LEVELS                                              QTPOS1B.403    
                DO I=TOP_ROW_START,TOP_ROW_START+ROW_LENGTH-1              QTPOS1B.404    
                  QTFIX(I,K)=NP_VALUE(K,2)                                 QTPOS1B.405    
                ENDDO                                                      QTPOS1B.406    
              ENDDO                                                        QTPOS1B.407    
            ENDIF                                                          QTPOS1B.408    
          ENDIF                                                            QTPOS1B.409    
                                                                           QTPOS1B.410    
! Set up arrays containing SP values                                       QTPOS1B.411    
          IF (at_base_of_LPG .AND. at_left_of_LPG) THEN                    QTPOS1B.412    
            DO K=1,Q_LEVELS                                                QTPOS1B.413    
              SP_VALUE(K,1)=QT(I_end_local,K)                              QTPOS1B.414    
              SP_VALUE(K,2)=QTFIX(I_end_local,K)                           QTPOS1B.415    
            ENDDO                                                          QTPOS1B.416    
          ELSE                                                             QTPOS1B.417    
            DO K=1,Q_LEVELS                                                QTPOS1B.418    
              SP_VALUE(K,1)=0.0                                            QTPOS1B.419    
              SP_VALUE(K,2)=0.0                                            QTPOS1B.420    
            ENDDO                                                          QTPOS1B.421    
          ENDIF                                                            QTPOS1B.422    
                                                                           QTPOS1B.423    
! Distribute the SP value over processors on polar row                     QTPOS1B.424    
          IF (at_base_of_LPG) THEN                                         QTPOS1B.425    
            CALL GCG_RSUM(Q_LEVELS,GC_ROW_GROUP,info,SP_VALUE(1,1))        QTPOS1B.426    
            DO K=1,Q_LEVELS                                                QTPOS1B.427    
              DO I=P_BOT_ROW_START,P_BOT_ROW_START+ROW_LENGTH-1            QTPOS1B.428    
                QT(I,K)=SP_VALUE(K,1)                                      QTPOS1B.429    
              ENDDO                                                        QTPOS1B.430    
            ENDDO                                                          QTPOS1B.431    
            IF (SF_QTFIX) THEN                                             QTPOS1B.432    
              CALL GCG_RSUM(Q_LEVELS,GC_ROW_GROUP,info,NP_VALUE(1,2))      QTPOS1B.433    
              DO K=1,Q_LEVELS                                              QTPOS1B.434    
                DO I=P_BOT_ROW_START,P_BOT_ROW_START+ROW_LENGTH-1          QTPOS1B.435    
                  QTFIX(I,K)=SP_VALUE(K,2)                                 QTPOS1B.436    
                ENDDO                                                      QTPOS1B.437    
              ENDDO                                                        QTPOS1B.438    
            ENDIF                                                          QTPOS1B.439    
          ENDIF                                                            QTPOS1B.440    
*ENDIF                                                                     QTPOS1B.441    
! Reset the pole rows of latitude                                          QTPOS1B.442    
                                                                           QTPOS1B.443    
*IF DEF,MPP                                                                QTPOS1B.444    
        IF (at_top_of_LPG .AND. at_right_of_LPG) THEN                      QTPOS1B.445    
*ENDIF                                                                     QTPOS1B.446    
          COS_P_LATITUDE(TOP_ROW_START+LAST_ROW_PT-1)=COS_P_LAT_NP         QTPOS1B.447    
          SEC_P_LATITUDE(TOP_ROW_START+LAST_ROW_PT-1)=SEC_P_LAT_NP         QTPOS1B.448    
*IF DEF,MPP                                                                QTPOS1B.449    
        ENDIF                                                              QTPOS1B.450    
                                                                           QTPOS1B.451    
        IF (at_base_of_LPG .AND. at_left_of_LPG) THEN                      QTPOS1B.452    
*ENDIF                                                                     QTPOS1B.453    
          COS_P_LATITUDE(P_BOT_ROW_START+FIRST_ROW_PT-1)=COS_P_LAT_SP      QTPOS1B.454    
          SEC_P_LATITUDE(P_BOT_ROW_START+FIRST_ROW_PT-1)=SEC_P_LAT_SP      QTPOS1B.455    
*IF DEF,MPP                                                                QTPOS1B.456    
        ENDIF                                                              QTPOS1B.457    
*ENDIF                                                                     QTPOS1B.458    
      ENDIF ! IF (.NOT. L_QT_POS_LOCAL)                                    QTPOS1B.459    
*ENDIF                                                                     QTPOS1B.460    
                                                                           QTPOS1B.461    
*IF DEF,T3E                                                                GPB0F405.198    
      CALL REMOVE_NEGATIVE_ZERO(QT,P_FIELD,Q_LEVELS)                       GPB0F405.199    
*ENDIF                                                                     GPB0F405.200    
*IF DEF,MPP                                                                QTPOS1B.462    
! Bring halos up to date                                                   QTPOS1B.463    
      CALL SWAPBOUNDS(QT,ROW_LENGTH,tot_P_ROWS,                            QTPOS1B.464    
     &                EW_Halo,NS_Halo,Q_LEVELS)                            QTPOS1B.465    
*ENDIF                                                                     QTPOS1B.466    
                                                                           QTPOS1B.467    
      RETURN                                                               QTPOS1B.468    
      END                                                                  QTPOS1B.469    
                                                                           QTPOS1B.470    
CLL   SUBROUTINE QT_POS -------------------------------------------        QTPOS1B.471    
CLL                                                                        QTPOS1B.472    
CLL   PURPOSE:   REMOVES NEGATIVE VALUES OF QT.  THERE ARE TWO ALTERNAT-   QTPOS1B.473    
CLL              IVE METHODS: METHOD 1 IS THE ORIGINAL SCHEME WHILE        QTPOS1B.474    
CLL              IF L_QT_POS_LOCAL = .TRUE. METHOD 2 IS USED.              QTPOS1B.475    
CLL              METHOD 2 IS DESIGNED TO ELIMINATE THE                     QTPOS1B.476    
CLL              VERY SLOW CLIMATE DRIFT IN QT IN UPPER MODEL LEVELS.      QTPOS1B.477    
CLL              IT IS UNNECESSARILY COMPLICATED (AND EXPENSIVE) FOR       QTPOS1B.478    
CLL              FORECAST USE.                                             QTPOS1B.479    
CLL                                                                        QTPOS1B.480    
CLL        METHOD 1:  ACCUMULATES TOTALS OF                                QTPOS1B.481    
CLL              OF NEGATIVE AND POSITIVE VALUES ON A LEVEL, ZEROING       QTPOS1B.482    
CLL              ALL NEGATIVE POINTS AND PROPORTIONALLY REMOVING THE       QTPOS1B.483    
CLL              SUM OF THE NEGATIVES FROM ALL POSITIVE POINTS.            QTPOS1B.484    
CLL        METHOD 2:  THE FOLLOWING METHOD IS APPLIED AT EACH LAYER:       QTPOS1B.485    
CLL            STEP 1 - LOOP ROUND ALL POINTS. IF QT IS NEGATIVE,          QTPOS1B.486    
CLL              (a):                                                      QTPOS1B.487    
CLL              SUM QT ROUND NEIGHBOURING POINTS WITHIN ONE ROW AND       QTPOS1B.488    
CLL              ONE COLUMN OF POINT IN CURRENT LAYER.  IF THIS VALUE      QTPOS1B.489    
CLL              IS POSITIVE AND LARGE ENOUGH,THE CENTRE NEGATIVE VALUE    QTPOS1B.490    
CLL              IS SET EQUAL TO ZERO AND THE NEIGHBOURS ARE SCALED        QTPOS1B.491    
CLL              TO CONSERVE MASS WEIGHTED QT. IF (a) DOES NOT WORK        QTPOS1B.492    
CLL              THEN (b):                                                 QTPOS1B.493    
CLL              EXTEND SUMMATION WITH A LARGER SEARCH RADIUS AND REPEAT   QTPOS1B.494    
CLL              PROCESS. (SEARCH RADIUS = 4 FOUND TO BE SUFFICIENT)       QTPOS1B.495    
CLL            STEP 2 - IF GLOBAL MODEL, CHECK FOR NEGATIVE QT AT POLES    QTPOS1B.496    
CLL              AND CORRECT USING ALL VALUES IN NEIGHBOURING ROW.         QTPOS1B.497    
CLL            STEP 3 - IF ANY NEGATIVE VALUES REMAIN (VERY UNCOMMON)      QTPOS1B.498    
CLL              PERFORM GLOBAL CORRECTION AS IN METHOD 1.                 QTPOS1B.499    
CLL                                                                        QTPOS1B.500    
CLL              IN BOTH METHODS IF THERE ARE INSUFFICIENT                 QTPOS1B.501    
CLL              POSITIVE VALUES WITHIN THE WHOLE LAYER AN ERROR           QTPOS1B.502    
CLL              MESSAGE IS PRODUCED.  PROGRAM CONTINUES IF L_NEG_Q        QTPOS1B.503    
CLL              IS TRUE IN WHICH CASE QT IS NOT CONSERVED.                QTPOS1B.504    
CLL   N.B.       DEFINITION OF RS SQUARED DELTAP MEANS THAT POSITIVE       QTPOS1B.505    
CLL              QT VALUES HAVE A NEGATIVE SIGN!                           QTPOS1B.506    
CLL                                                                        QTPOS1B.507    
CLL                                                                        QTPOS1B.508    
CLL   NOT SUITABLE FOR SINGLE COLUMN USE.                                  QTPOS1B.509    
CLL   VERSION FOR CRAY Y-MP                                                QTPOS1B.510    
CLL                                                                        QTPOS1B.511    
CLL  MM, SB     <- PROGRAMMER OF SOME OR ALL OF PREVIOUS CODE OR CHANGES   QTPOS1B.512    
CLL                                                                        QTPOS1B.513    
CLL  MODEL            MODIFICATION HISTORY FROM MODEL VERSION 4.2:         ATJ0F403.225    
CLL VERSION  DATE                                                          QTPOS1B.515    
!LL   4.3   14/04/97  Minor updates in line with QTPOS1A.  T Johns         ATJ0F403.226    
CLL                                                                        QTPOS1B.516    
CLL   PROGRAMMING STANDARD: UNIFIED MODEL DOCUMENTATION PAPER NO. 4,       QTPOS1B.517    
CLL                         STANDARD A. VERSION 2, DATED 18/01/90          QTPOS1B.518    
CLL                                                                        QTPOS1B.519    
CLL   SYSTEM COMPONENTS COVERED: P 191                                     QTPOS1B.520    
CLL                                                                        QTPOS1B.521    
CLL   SYSTEM TASK: P1                                                      QTPOS1B.522    
CLL                                                                        QTPOS1B.523    
CLL   DOCUMENTATION:       SECTION 3.7                                     QTPOS1B.524    
CLL                        IN UNIFIED MODEL DOCUMENTATION PAPER            QTPOS1B.525    
CLL                        NO. 10 M.J.P.CULLEN,T.DAVIES AND M.H.MAWSON     QTPOS1B.526    
CLL                        VERSION 10, DATED 10/09/90.                     QTPOS1B.527    
CLLEND-------------------------------------------------------------        QTPOS1B.528    
                                                                           QTPOS1B.529    
C*L   ARGUMENTS:---------------------------------------------------        QTPOS1B.530    

      SUBROUTINE QT_POS                                                     4QTPOS1B.531    
     1                 (QT,ROW_LENGTH,P_FIELD,ISTART,IEND,Q_LEVELS,        QTPOS1B.532    
     2                  ERROR_CODE,ERROR_MESSAGE,FACTOR,FAILURE,           QTPOS1B.533    
     3                  L_NEG_QT,L_QT_POS_LOCAL)                           QTPOS1B.534    
                                                                           QTPOS1B.535    
      IMPLICIT NONE                                                        QTPOS1B.536    
                                                                           QTPOS1B.537    
      INTEGER                                                              QTPOS1B.538    
     &  ROW_LENGTH         ! IN : number of points per row                 QTPOS1B.539    
     &, P_FIELD            ! IN : size of horizontal P field               QTPOS1B.540    
     &, ISTART             ! IN : first point to process                   QTPOS1B.541    
     &, IEND               ! IN : last point to process                    QTPOS1B.542    
     &, Q_LEVELS           ! IN : number of (moist) levels to process      QTPOS1B.543    
                                                                           QTPOS1B.544    
      LOGICAL                                                              QTPOS1B.545    
     &  L_NEG_QT           ! IN : switch to continue in event of           QTPOS1B.546    
!                          !      excessive negative QTs (non-conserv.)    QTPOS1B.547    
     &, L_QT_POS_LOCAL     ! IN : perform local corrections (method 2)     QTPOS1B.548    
                                                                           QTPOS1B.549    
      REAL                                                                 QTPOS1B.550    
     &  QT(P_FIELD,Q_LEVELS) ! IN/OUT : QT field to update                 QTPOS1B.551    
                                                                           QTPOS1B.552    
      INTEGER                                                              QTPOS1B.553    
     &  ERROR_CODE         ! OUT : 1 on exit if error detected             QTPOS1B.554    
                                                                           QTPOS1B.555    
      CHARACTER*(80)                                                       QTPOS1B.556    
     &  ERROR_MESSAGE      ! OUT : Error message if ERROR_CODE not zero    QTPOS1B.557    
                                                                           QTPOS1B.558    
      REAL                                                                 QTPOS1B.559    
     &  FACTOR(Q_LEVELS)   ! OUT : Scaling factor for QT                   QTPOS1B.560    
                                                                           QTPOS1B.561    
      LOGICAL                                                              QTPOS1B.562    
     &  FAILURE(Q_LEVELS)  ! OUT : Indicates a level needs rescaling       QTPOS1B.563    
                                                                           QTPOS1B.564    
C*---------------------------------------------------------------------    QTPOS1B.565    
                                                                           QTPOS1B.566    
C*L   DEFINE ARRAYS AND VARIABLES USED IN THIS ROUTINE-----------------    QTPOS1B.567    
C DEFINE LOCAL ARRAYS: NONE ARE REQUIRED                                   QTPOS1B.568    
C*---------------------------------------------------------------------    QTPOS1B.569    
C DEFINE LOCAL VARIABLES                                                   QTPOS1B.570    
                                                                           QTPOS1B.571    
      INTEGER                                                              QTPOS1B.572    
     *  I,K,II,JJ,NN,J   ! Loop variables                                  QTPOS1B.573    
     *, IN,IS,JE,JW        ! Row/column indices to north/south/east/west   QTPOS1B.574    
     *, POINTER            ! Array pointer                                 QTPOS1B.575    
     *, ROW                ! Row number counting from 0                    QTPOS1B.576    
     *, NPT                ! Point in row counting from 1                  QTPOS1B.577    
     *, ROWS               ! Number of rows                                QTPOS1B.578    
     *, MAX_SEARCH         ! Maximum search radius (points)                QTPOS1B.579    
     *, NO_NEG             ! Count of number negative in step 3 & 1        QTPOS1B.580    
     *, NO_NEG2            ! Count of number negative                      QTPOS1B.581    
     *, NINDEX(P_FIELD)    ! locations of negative q                       QTPOS1B.582    
      REAL                                                                 QTPOS1B.583    
     *  SUM_NEIGHBOURS     ! Sum of QT around neighbours                   QTPOS1B.584    
     *, SUM_POSITIVE       ! Global sum of positive QT                     QTPOS1B.585    
     *, SUM_NEGATIVE       ! Global sum of negative QT                     QTPOS1B.586    
     *, TEMP1,TEMP2        ! Temporary sums                                QTPOS1B.587    
                                                                           QTPOS1B.588    
C*L   NO EXTERNAL SUBROUTINE CALLS:------------------------------------    QTPOS1B.589    
C*---------------------------------------------------------------------    QTPOS1B.590    
CL    MAXIMUM VECTOR LENGTH ASSUMED IS P_FIELD.                            QTPOS1B.591    
CL---------------------------------------------------------------------    QTPOS1B.592    
CL    INTERNAL STRUCTURE.                                                  QTPOS1B.593    
CL---------------------------------------------------------------------    QTPOS1B.594    
CL                                                                         QTPOS1B.595    
      ROWS=P_FIELD/ROW_LENGTH                                              QTPOS1B.596    
      MAX_SEARCH=4                                                         QTPOS1B.597    
                                                                           QTPOS1B.598    
CL LOOP OVER LEVELS                                                        QTPOS1B.599    
      DO K= 1,Q_LEVELS                                                     QTPOS1B.600    
        FAILURE(K)=.FALSE.                                                 QTPOS1B.601    
        FACTOR(K)=0.0                                                      QTPOS1B.602    
        IF (L_QT_POS_LOCAL) THEN                                           QTPOS1B.603    
CL---------------------------------------------------------------------    QTPOS1B.604    
CL    METHOD 2                                                             QTPOS1B.605    
CL---------------------------------------------------------------------    QTPOS1B.606    
          NO_NEG=0                                                         QTPOS1B.607    
                                                                           QTPOS1B.608    
*IF DEF,GLOBAL                                                             QTPOS1B.609    
C   Loop round non-polar points                                            QTPOS1B.610    
          DO I=ROW_LENGTH+1,P_FIELD-ROW_LENGTH                             QTPOS1B.611    
*ELSE                                                                      QTPOS1B.612    
          DO I=1,P_FIELD                                                   QTPOS1B.613    
*ENDIF                                                                     QTPOS1B.614    
            IF(QT(I,K).GT.0.0) THEN                                        QTPOS1B.615    
              NO_NEG=NO_NEG+1                                              QTPOS1B.616    
              NINDEX(NO_NEG)=I                                             QTPOS1B.617    
            ENDIF                                                          QTPOS1B.618    
          ENDDO                                                            QTPOS1B.619    
                                                                           QTPOS1B.620    
          DO J=1,NO_NEG                                                    QTPOS1B.621    
            I=NINDEX(J)                                                    QTPOS1B.622    
CL                                                                         QTPOS1B.623    
CL---------------------------------------------------------------------    QTPOS1B.624    
CL    SECTION 2.    STEP 1. WHERE QT IS NEGATIVE LOOP ROUND                QTPOS1B.625    
CL              8 NEIGHBOURS TO PERFORM LOCAL CORRECTION.  THIS OCCURS     QTPOS1B.626    
CL              AT APPROXIMATELY 2% OF THE POINTS.  THE FOLLOWING CODE     QTPOS1B.627    
CL              IN SECTIONS 2 & 3 IS RATHER LONG-WINDED, BUT CPU IS        QTPOS1B.628    
CL              SIGNIFICANTLY REDUCED OVER A SIMPLER APPROACH.             QTPOS1B.629    
CL---------------------------------------------------------------------    QTPOS1B.630    
C   If QT is negative find 8 neighbours ...                                QTPOS1B.631    
              ROW=(I-1)/ROW_LENGTH                                         QTPOS1B.632    
              NPT=I-ROW*ROW_LENGTH                                         QTPOS1B.633    
              JE=1                                                         QTPOS1B.634    
              JW=-1                                                        QTPOS1B.635    
*IF DEF,GLOBAL                                                             QTPOS1B.636    
              IF(NPT.EQ.ROW_LENGTH) JE=1-ROW_LENGTH                        QTPOS1B.637    
              IF(NPT.EQ.1) JW=ROW_LENGTH-1                                 QTPOS1B.638    
*ELSE                                                                      QTPOS1B.639    
C   Abandon search for neighbours on LAM boundary                          QTPOS1B.640    
              IF(NPT.EQ.ROW_LENGTH.OR.NPT.EQ.1) THEN                       QTPOS1B.641    
                FAILURE(K)=.TRUE.                                          QTPOS1B.642    
                GOTO 200                                                   QTPOS1B.643    
              END IF                                                       QTPOS1B.644    
*ENDIF                                                                     QTPOS1B.645    
C   ... and sum neighbouring 8 values                                      QTPOS1B.646    
              SUM_NEIGHBOURS = QT(I+JE,K) + QT(I+JW,K)                     QTPOS1B.647    
              IF(ROW.GT.0) THEN                                            QTPOS1B.648    
                II=I-ROW_LENGTH                                            QTPOS1B.649    
                SUM_NEIGHBOURS = SUM_NEIGHBOURS +                          QTPOS1B.650    
     *                           QT(II,K) + QT(II+JE,K) + QT(II+JW,K)      QTPOS1B.651    
              END IF                                                       QTPOS1B.652    
              IF(ROW.LT.ROWS-1) THEN                                       QTPOS1B.653    
                II=I+ROW_LENGTH                                            QTPOS1B.654    
                SUM_NEIGHBOURS = SUM_NEIGHBOURS +                          QTPOS1B.655    
     *                           QT(II,K) + QT(II+JE,K) + QT(II+JW,K)      QTPOS1B.656    
              END IF                                                       QTPOS1B.657    
C                                                                          QTPOS1B.658    
C   If possible set QT=0 at point with negative value and adjust           QTPOS1B.659    
C   neighbouring values by FACTOR to conserve QT                           QTPOS1B.660    
              IF(SUM_NEIGHBOURS.LT.0.0) THEN                               QTPOS1B.661    
                FACTOR(K)=1.0+QT(I,K)/SUM_NEIGHBOURS                       QTPOS1B.662    
                IF(FACTOR(K).GE.0.0) THEN                                  QTPOS1B.663    
                  QT(I,K)=0.0                                              QTPOS1B.664    
                  QT(I+JE,K) = QT(I+JE,K)*FACTOR(K)                        QTPOS1B.665    
                  QT(I+JW,K) = QT(I+JW,K)*FACTOR(K)                        QTPOS1B.666    
                  IF(ROW.GT.0) THEN                                        QTPOS1B.667    
                    II=I-ROW_LENGTH                                        QTPOS1B.668    
                    QT(II,K) = QT(II,K)*FACTOR(K)                          QTPOS1B.669    
                    QT(II+JE,K) = QT(II+JE,K)*FACTOR(K)                    QTPOS1B.670    
                    QT(II+JW,K) = QT(II+JW,K)*FACTOR(K)                    QTPOS1B.671    
                  END IF                                                   QTPOS1B.672    
                  IF(ROW.LT.ROWS-1) THEN                                   QTPOS1B.673    
                    II=I+ROW_LENGTH                                        QTPOS1B.674    
                    QT(II,K) = QT(II,K)*FACTOR(K)                          QTPOS1B.675    
                    QT(II+JE,K) = QT(II+JE,K)*FACTOR(K)                    QTPOS1B.676    
                    QT(II+JW,K) = QT(II+JW,K)*FACTOR(K)                    QTPOS1B.677    
                  END IF                                                   QTPOS1B.678    
                END IF                                                     QTPOS1B.679    
              END IF                                                       QTPOS1B.680    
*IF -DEF,GLOBAL                                                            QTPOS1B.681    
  200       CONTINUE       ! LAM only                                      QTPOS1B.682    
*ENDIF                                                                     QTPOS1B.683    
          END DO       ! loop over negative grid points                    QTPOS1B.684    
                                                                           QTPOS1B.685    
! recalculate number of negative gridpoints                                QTPOS1B.686    
          NO_NEG2=0                                                        QTPOS1B.687    
          DO J=1,NO_NEG                                                    QTPOS1B.688    
            I=NINDEX(J)                                                    QTPOS1B.689    
            IF(QT(I,K).GT.0.0) THEN                                        QTPOS1B.690    
              NO_NEG2=NO_NEG2+1                                            QTPOS1B.691    
              NINDEX(NO_NEG2)=I                                            QTPOS1B.692    
            ENDIF                                                          QTPOS1B.693    
          ENDDO                                                            QTPOS1B.694    
!         WRITE(6,*),' Step 1 level ',k,' no neg ',no_NEG                  ATJ0F403.227    
!         WRITE(6,*),' Step 2 level ',k,' no neg ',no_NEG2                 ATJ0F403.228    
CL                                                                         QTPOS1B.697    
CL---------------------------------------------------------------------    QTPOS1B.698    
CL    SECTION 3.    STEP 2. WHERE QT IS NEGATIVE AND STEP 1 CORRECTIONS    QTPOS1B.699    
CL              HAVE FAILED, ATTEMPT LOCAL CORRECTION BY EXTENDING         QTPOS1B.700    
CL              SEARCH TO 24, 48 OR 80 SURROUNDING POINTS (SEARCH RADIUS   QTPOS1B.701    
CL              EQUALS 2, 3 OR 4), OR FURTHER.  THIS OCCURS ABOUT          QTPOS1B.702    
CL              0.1% OF THE TIME AT CLIMATE RESOLUTION AND MUCH MORE       QTPOS1B.703    
CL              FREQUENTLY AT FORECAST RESOLUTION.  ALTHOUGH SMALL IN      QTPOS1B.704    
CL              NUMBER, LOCAL CORRECTIONS AT THESE POINTS ARE REQUIRED     QTPOS1B.705    
CL              TO SLOW DOWN THE LOSS OF STRATOSPHERIC MOISTURE THAT       QTPOS1B.706    
CL              OCCURS ON CLIMATE TIME SCALES.                             QTPOS1B.707    
CL---------------------------------------------------------------------    QTPOS1B.708    
         IF (NO_NEG2.GT.0) THEN         ! only do this if required         QTPOS1B.709    
           DO J=1,NO_NEG2                                                  QTPOS1B.710    
            I=NINDEX(J)                                                    QTPOS1B.711    
              ROW=(I-1)/ROW_LENGTH                                         QTPOS1B.712    
              NPT=I-ROW*ROW_LENGTH                                         QTPOS1B.713    
*IF -DEF,GLOBAL                                                            QTPOS1B.714    
C   Abandon search for neighbours on LAM boundary                          QTPOS1B.715    
              IF(NPT.EQ.ROW_LENGTH.OR.NPT.EQ.1) THEN                       QTPOS1B.716    
                FAILURE(K)=.TRUE.                                          QTPOS1B.717    
                GOTO 300                                                   QTPOS1B.718    
              END IF                                                       QTPOS1B.719    
*ENDIF                                                                     QTPOS1B.720    
              DO NN=2,MAX_SEARCH                                           QTPOS1B.721    
                JW=NPT-NN                                                  QTPOS1B.722    
                JE=NPT+NN                                                  QTPOS1B.723    
*IF -DEF,GLOBAL                                                            QTPOS1B.724    
                IF(JE.GT.ROW_LENGTH) JE=ROW_LENGTH                         QTPOS1B.725    
                IF(JW.LT.1) JW=1                                           QTPOS1B.726    
*ENDIF                                                                     QTPOS1B.727    
                IN=ROW-NN                                                  QTPOS1B.728    
                IS=ROW+NN                                                  QTPOS1B.729    
                IF(IN.LT.0) IN=0                                           QTPOS1B.730    
                IF(IS.GE.ROWS) IS=ROWS-1                                   QTPOS1B.731    
C   ... and sum QT at neighbouring points                                  QTPOS1B.732    
                SUM_NEIGHBOURS=-QT(I,K)                                    QTPOS1B.733    
                DO II=IN,IS                                                QTPOS1B.734    
                  DO JJ=JW,JE                                              QTPOS1B.735    
                    POINTER = II*ROW_LENGTH + JJ                           QTPOS1B.736    
*IF DEF,GLOBAL                                                             QTPOS1B.737    
                    IF(JJ.GT.ROW_LENGTH) POINTER = POINTER-ROW_LENGTH      QTPOS1B.738    
                    IF(JJ.LT.1) POINTER = POINTER+ROW_LENGTH               QTPOS1B.739    
*ENDIF                                                                     QTPOS1B.740    
                    SUM_NEIGHBOURS = SUM_NEIGHBOURS + QT(POINTER,K)        QTPOS1B.741    
                  END DO                                                   QTPOS1B.742    
                END DO                                                     QTPOS1B.743    
C                                                                          QTPOS1B.744    
C   If possible set QT=0 at point with negative value and adjust           QTPOS1B.745    
C   neighbouring values by FACTOR to conserve QT                           QTPOS1B.746    
                IF(SUM_NEIGHBOURS.LT.0.0) THEN                             QTPOS1B.747    
                  FACTOR(K)=1.0+QT(I,K)/SUM_NEIGHBOURS                     QTPOS1B.748    
                  IF(FACTOR(K).GE.0.0) THEN                                QTPOS1B.749    
                    QT(I,K)=0.0                                            QTPOS1B.750    
                    DO II=IN,IS                                            QTPOS1B.751    
                      DO JJ=JW,JE                                          QTPOS1B.752    
                        POINTER = II*ROW_LENGTH + JJ                       QTPOS1B.753    
*IF DEF,GLOBAL                                                             QTPOS1B.754    
                        IF(JJ.GT.ROW_LENGTH) POINTER =POINTER-ROW_LENGTH   QTPOS1B.755    
                        IF(JJ.LT.1) POINTER = POINTER+ROW_LENGTH           QTPOS1B.756    
*ENDIF                                                                     QTPOS1B.757    
                        QT(POINTER,K) = QT(POINTER,K)*FACTOR(K)            QTPOS1B.758    
                      END DO                                               QTPOS1B.759    
                    END DO                                                 QTPOS1B.760    
                    GOTO 300    ! Successfull correction performed         QTPOS1B.761    
                  END IF                                                   QTPOS1B.762    
                END IF                                                     QTPOS1B.763    
              END DO     ! End loop on search radius (NN)                  QTPOS1B.764    
C                                                                          QTPOS1B.765    
C  No correction possible within maximum search radius (=4)                QTPOS1B.766    
              FAILURE(K)=.TRUE.                                            QTPOS1B.767    
  300       CONTINUE                                                       QTPOS1B.768    
          END DO         ! End loop on points (I)                          QTPOS1B.769    
         END IF          ! end step 2                                      QTPOS1B.770    
                                                                           QTPOS1B.771    
*IF DEF,GLOBAL                                                             QTPOS1B.772    
C                                                                          QTPOS1B.773    
C  Loop round polar points                                                 QTPOS1B.774    
CL                                                                         QTPOS1B.775    
CL---------------------------------------------------------------------    QTPOS1B.776    
CL    SECTION 4.     STEP 2. IF GLOBAL MODEL DEAL WITH POLES IN            QTPOS1B.777    
CL                   SIMILAR FASHION EXCEPT WORKING ROW AT A TIME.         QTPOS1B.778    
CL---------------------------------------------------------------------    QTPOS1B.779    
C  Set all ploar points equal                                              QTPOS1B.780    
          TEMP1=0.0                                                        QTPOS1B.781    
          TEMP2=0.0                                                        QTPOS1B.782    
          DO I=1,ROW_LENGTH                                                QTPOS1B.783    
            TEMP1=TEMP1+QT(I,K)                                            QTPOS1B.784    
            TEMP2=TEMP2+QT(P_FIELD-ROW_LENGTH+I,K)                         QTPOS1B.785    
          END DO                                                           QTPOS1B.786    
          DO I=1,ROW_LENGTH                                                QTPOS1B.787    
            QT(I,K)=TEMP1/ROW_LENGTH                                       QTPOS1B.788    
            QT(P_FIELD-ROW_LENGTH+I,K)=TEMP2/ROW_LENGTH                    QTPOS1B.789    
          END DO                                                           QTPOS1B.790    
C  Check for negative QT at north pole                                     QTPOS1B.791    
          IF(QT(1,K).GT.0.0) THEN                                          QTPOS1B.792    
            SUM_NEIGHBOURS=0.0                                             QTPOS1B.793    
            DO I=1,ROW_LENGTH                                              QTPOS1B.794    
              SUM_NEIGHBOURS = SUM_NEIGHBOURS + QT(ROW_LENGTH+I,K)         QTPOS1B.795    
            END DO                                                         QTPOS1B.796    
            IF(SUM_NEIGHBOURS.LT.0.0) THEN                                 QTPOS1B.797    
              FACTOR(K)=1.0+QT(1,K)*ROW_LENGTH/SUM_NEIGHBOURS              QTPOS1B.798    
              IF(FACTOR(K).GE.0.0) THEN                                    QTPOS1B.799    
                DO I=1,ROW_LENGTH                                          QTPOS1B.800    
                  QT(I,K)=0.0                                              QTPOS1B.801    
                  QT(ROW_LENGTH+I,K)=QT(ROW_LENGTH+I,K)*FACTOR(K)          QTPOS1B.802    
                END DO                                                     QTPOS1B.803    
                GOTO 400    ! Correction performed                         QTPOS1B.804    
              END IF                                                       QTPOS1B.805    
            END IF                                                         QTPOS1B.806    
C                                                                          QTPOS1B.807    
C  No correction possible: global QT in layer is negative                  QTPOS1B.808    
            FAILURE(K)=.TRUE.                                              QTPOS1B.809    
                                                                           QTPOS1B.810    
          END IF                                                           QTPOS1B.811    
  400     CONTINUE                                                         QTPOS1B.812    
C  Check for negative QT at south pole                                     QTPOS1B.813    
          IF(QT(P_FIELD-1,K).GT.0.0) THEN                                  QTPOS1B.814    
            SUM_NEIGHBOURS=0.0                                             QTPOS1B.815    
            DO I=1,ROW_LENGTH                                              QTPOS1B.816    
              SUM_NEIGHBOURS=SUM_NEIGHBOURS+QT(P_FIELD+1-ROW_LENGTH-I,K)   QTPOS1B.817    
            END DO                                                         QTPOS1B.818    
            IF(SUM_NEIGHBOURS.LT.0.0) THEN                                 QTPOS1B.819    
              FACTOR(K)=1.0+QT(P_FIELD-1,K)*ROW_LENGTH/SUM_NEIGHBOURS      QTPOS1B.820    
              IF(FACTOR(K).GE.0.0) THEN                                    QTPOS1B.821    
CDIR$ IVDEP                                                                QTPOS1B.822    
                DO I=1,ROW_LENGTH                                          QTPOS1B.823    
                  QT(P_FIELD-ROW_LENGTH+I,K)=0.0                           QTPOS1B.824    
                  QT(P_FIELD+1-ROW_LENGTH-I,K) =                           QTPOS1B.825    
     *                       QT(P_FIELD+1-ROW_LENGTH-I,K)*FACTOR(K)        QTPOS1B.826    
                END DO                                                     QTPOS1B.827    
                GOTO 500    ! Correction performed                         QTPOS1B.828    
              END IF                                                       QTPOS1B.829    
            END IF                                                         QTPOS1B.830    
C                                                                          QTPOS1B.831    
C  No correction possible: global QT in layer is negative                  QTPOS1B.832    
            FAILURE(K)=.TRUE.                                              QTPOS1B.833    
          END IF                                                           QTPOS1B.834    
  500     CONTINUE                                                         QTPOS1B.835    
*ENDIF                                                                     QTPOS1B.836    
        END IF          !  End main section for method 1                   QTPOS1B.837    
CL---------------------------------------------------------------------    QTPOS1B.838    
CL    SECTION 5.     METHOD 1.  ALSO METHOD 2 STEP 3.                      QTPOS1B.839    
CL              REMOVE REMAINING NEGATIVE QT BY SUMMING                    QTPOS1B.840    
CL              NEGATIVE AND POSITIVE VALUES SEPARATELY ON THE             QTPOS1B.841    
CL              LEVEL AND PERFORMING A GLOBAL CORRECTION.                  QTPOS1B.842    
CL---------------------------------------------------------------------    QTPOS1B.843    
        IF(FAILURE(K).OR..NOT.L_QT_POS_LOCAL) THEN                         QTPOS1B.844    
          NO_NEG=0            !zero count of -ve points                    QTPOS1B.845    
          SUM_POSITIVE = 0.0                                               QTPOS1B.846    
          SUM_NEGATIVE = 0.0                                               QTPOS1B.847    
          DO I=ISTART,IEND                                                 QTPOS1B.848    
            IF(QT(I,K).GT.0.)THEN                                          QTPOS1B.849    
              NO_NEG=NO_NEG+1           !count number of -ve points        QTPOS1B.850    
              SUM_NEGATIVE = SUM_NEGATIVE + QT(I,K)                        QTPOS1B.851    
              QT(I,K) = 0.0                                                QTPOS1B.852    
            ELSE                                                           QTPOS1B.853    
              SUM_POSITIVE = SUM_POSITIVE + QT(I,K)                        QTPOS1B.854    
            END IF                                                         QTPOS1B.855    
          END DO                                                           QTPOS1B.856    
!         WRITE(6,*),' GLOBAL correction level ',k,' no neg ',no_NEG       ATJ0F403.229    
                                                                           QTPOS1B.858    
CL If a negative value is found re-scale all positive points               QTPOS1B.859    
CL                                                                         QTPOS1B.860    
          IF (SUM_NEGATIVE.GT.0.0) THEN                                    QTPOS1B.861    
            FAILURE(K)=.TRUE.                                              QTPOS1B.862    
            FACTOR(K) = 1. + SUM_NEGATIVE/SUM_POSITIVE                     QTPOS1B.863    
            IF(L_NEG_QT.AND.FACTOR(K).LT.0)THEN                            QTPOS1B.864    
CL Skip abort if L_NEG_QT is true                                          QTPOS1B.865    
              WRITE(6,*)' QT_POS : Mass weighted QT summed over level ',   QTPOS1B.866    
     *            K,' was negative.  WARNING: QT not conserved'            QTPOS1B.867    
              WRITE(6,*) NO_NEG,' points were -ve and the scaling ',       QTPOS1B.868    
     *           'factor has been reset to 1'                              QTPOS1B.869    
              FACTOR(K) = 1.0                                              QTPOS1B.870    
            ENDIF                                                          QTPOS1B.871    
            IF(FACTOR(K).LT.0.0) THEN                                      QTPOS1B.872    
              ERROR_CODE = 1                                               QTPOS1B.873    
              ERROR_MESSAGE=                                               QTPOS1B.874    
     *  'QT_POS   : MASS-WEIGHTED QT SUMMED OVER A LEVEL WAS NEGATIVE.'    QTPOS1B.875    
              RETURN                                                       QTPOS1B.876    
            END IF                                                         QTPOS1B.877    
          ELSE                                                             QTPOS1B.878    
            FAILURE(K)=.FALSE.                                             QTPOS1B.879    
          END IF                                                           QTPOS1B.880    
        END IF                                                             QTPOS1B.881    
                                                                           QTPOS1B.882    
CL END LOOP OVER LEVELS                                                    QTPOS1B.883    
      END DO                                                               QTPOS1B.884    
                                                                           QTPOS1B.885    
                                                                           QTPOS1B.886    
CL    END OF ROUTINE QT_POS                                                QTPOS1B.887    
                                                                           QTPOS1B.888    
      RETURN                                                               QTPOS1B.889    
      END                                                                  QTPOS1B.890    
*ENDIF                                                                     QTPOS1B.891