C *****************************COPYRIGHT******************************             RAD_ALL.25250 
C (c) CROWN COPYRIGHT 1996, METEOROLOGICAL OFFICE, All Rights Reserved.            RAD_ALL.25251 
C                                                                                  RAD_ALL.25252 
C Use, duplication or disclosure of this code is subject to the                    RAD_ALL.25253 
C restrictions as set forth in the contract.                                       RAD_ALL.25254 
C                                                                                  RAD_ALL.25255 
C                Meteorological Office                                             RAD_ALL.25256 
C                London Road                                                       RAD_ALL.25257 
C                BRACKNELL                                                         RAD_ALL.25258 
C                Berkshire UK                                                      RAD_ALL.25259 
C                RG12 2SZ                                                          RAD_ALL.25260 
C                                                                                  RAD_ALL.25261 
C If no contract has been raised with this copy of the code, the use,              RAD_ALL.25262 
C duplication or disclosure of it is strictly prohibited.  Permission              RAD_ALL.25263 
C to do so must first be obtained in writing from the Head of Numerical            RAD_ALL.25264 
C Modelling at the above address.                                                  RAD_ALL.25265 
C ******************************COPYRIGHT******************************            RAD_ALL.25266 
CLL  SUBROUTINE TROPIN------------------------------------------------             RAD_ALL.25267 
CLL                                                                                RAD_ALL.25268 
CLL  Purpose:  Finds the tropopause & returns index of where it is.                RAD_ALL.25269 
CLL      Not suitable for single-column model use, as it does                      RAD_ALL.25270 
CLL             horizontal filling-in, & includes dynamical allocation.            RAD_ALL.25271 
CLL                                                                                RAD_ALL.25272 
CLL    Based on routine TROP, but                                                  RAD_ALL.25273 
CLL      1) taking temperature rather than theta as input                          RAD_ALL.25274 
CLL      2) rather than returning the temperature, pressure and height             RAD_ALL.25275 
CLL      of a continuously-varying tropopause found by extrapolating               RAD_ALL.25276 
CLL      lapse rates, it returns the index of the adjacent layer                   RAD_ALL.25277 
CLL      boundary (N meaning the bottom of the Nth model layer)                    RAD_ALL.25278 
CLL      3) "filling in" where no tropopause is found.                             RAD_ALL.25279 
CLL                                                                                RAD_ALL.25280 
CLL        Author:  William Ingram                                                 RAD_ALL.25281 
CLL                                                                                RAD_ALL.25282 
CLL  Model                Modification history:                                    RAD_ALL.25283 
CLL version  Date                                                                  RAD_ALL.25284 
CLL   4.2   23/9/96       New code, based on TROP                                  RAD_ALL.25285 
CLL                                                                                RAD_ALL.25286 
CLL    Note that the definition of "tropopause" matches TROP, which                RAD_ALL.25287 
CLL      is not quite the WMO definition, though the same critical                 RAD_ALL.25288 
CLL      lapse rate is used (unless comdeck C_LAPSE is altered).                   RAD_ALL.25289 
CLL      For details of the interpolation assumptions, see UMDP S1                 RAD_ALL.25290 
CLL      section 3.2.2 or Swinbank & Wilson (1990: SRFRTN 48).                     RAD_ALL.25291 
CLL      Any physical changes to one routine should be considered for              RAD_ALL.25292 
CLL                                         mirroring in the other.                RAD_ALL.25293 
CLLEND----------------------------------------------------------------             RAD_ALL.25294 
C                                                                                  RAD_ALL.25295 
C*L  Arguments:-------------------------------------------------------             RAD_ALL.25296 
      SUBROUTINE TROPIN(PSTAR, T, P_EXNER_HALF, IT, L1, L2, ROW_LENGTH,            RAD_ALL.25297 
     &     P_LEVELS, MIN_TROP_LEVEL, MAX_TROP_LEVEL, AKH, BKH, WRAP)               RAD_ALL.25298 
                                                                                   RAD_ALL.25299 
      IMPLICIT NONE                                                                RAD_ALL.25300 
                                                                                   RAD_ALL.25301 
      INTEGER!, INTENT (IN)                                                        RAD_ALL.25302 
     &     L1,                          !   Number of points in arrays             RAD_ALL.25303 
     &     L2,                          !   Number of points to process            RAD_ALL.25304 
     &     ROW_LENGTH,                  !   Number of points per row               RAD_ALL.25305 
     &     P_LEVELS,                    !   Number of model levels                 RAD_ALL.25306 
     &     MIN_TROP_LEVEL,                                                         RAD_ALL.25307 
     &     MAX_TROP_LEVEL                                                          RAD_ALL.25308 
C     ! Limits on where the tropopause may be deemed to be - between               RAD_ALL.25309 
C     !  the MIN_TROP_LEVELth and MAX_TROP_LEVELth layers (with the                RAD_ALL.25310 
C     !  convention used here for layer boundaries, the actual index               RAD_ALL.25311 
C     !  IT returned has MIN_TROP_LEVEL < IT =< MAX_TROP_LEVEL.)                   RAD_ALL.25312 
                                                                                   RAD_ALL.25313 
      REAL!, INTENT(IN)                                                            RAD_ALL.25314 
     &     PSTAR(L1),                   !  Surface pressure                        RAD_ALL.25315 
     &     T(L1,P_LEVELS),              !  Temperature at layer centres            RAD_ALL.25316 
     &     P_EXNER_HALF(L1,P_LEVELS+1), !  Pexner at layer boundaries              RAD_ALL.25317 
     &     AKH(P_LEVELS+1),             !  Hybrid coordinate A & B                 RAD_ALL.25318 
     &     BKH(P_LEVELS+1)              !   values for layer boundaries            RAD_ALL.25319 
                                                                                   RAD_ALL.25320 
      LOGICAL!, INTENT (IN)                                                        RAD_ALL.25321 
     &     WRAP                                                                    RAD_ALL.25322 
C     ! Do the rows wrap round (so that the last point of a row is                 RAD_ALL.25323 
C     !  geographically beside the first point of the same row) ?                  RAD_ALL.25324 
                                                                                   RAD_ALL.25325 
      INTEGER!, INTENT (OUT)                                                       RAD_ALL.25326 
     &     IT(L1)                                                                  RAD_ALL.25327 
C     ! Integer indexing the tropopause, taken to be @ a layer boundary            RAD_ALL.25328 
C     !   with the convention that N means the bottom of layer N.                  RAD_ALL.25329 
                                                                                   RAD_ALL.25330 
C Workspace usage:-----------------------------------------------------            RAD_ALL.25331 
                                                                                   RAD_ALL.25332 
      REAL LAPSE_RATE(L2,                                                          RAD_ALL.25333 
     & MIN_TROP_LEVEL+1:MAX_TROP_LEVEL+1-1/(1+P_LEVELS-MAX_TROP_LEVEL))            RAD_ALL.25334 
C     ! Lapse rate between layer centres (with top level MAX_TROP_LEVEL            RAD_ALL.25335 
C     !    if MAX_TROP_LEVEL = P_LEVELS, MAX_TROP_LEVEL+1 otherwise,               RAD_ALL.25336 
C     !    assuming MAX_TROP_LEVEL =< P_LEVELS as it should be).                   RAD_ALL.25337 
                                                                                   RAD_ALL.25338 
      LOGICAL LTROP(L2)                                                            RAD_ALL.25339 
C     !  Logical array which indicates whether we are still seeking a              RAD_ALL.25340 
C     !    tropopause (if not, it has already been found lower down)               RAD_ALL.25341 
                                                                                   RAD_ALL.25342 
C*---------------------------------------------------------------------            RAD_ALL.25343 
C Define local variables:----------------------------------------------            RAD_ALL.25344 
                                                                                   RAD_ALL.25345 
      INTEGER I, J,         ! Loopers over level & point                           RAD_ALL.25346 
     &     POINT,           ! Point counter at ends of rows                        RAD_ALL.25347 
     &     JP1,                                                                    RAD_ALL.25348 
C     !  J+1, except where this would cause out-of-bounds reference                RAD_ALL.25349 
     &     NNEIGH,          ! Number of well-defined tropopauses among             RAD_ALL.25350 
C     ! the 8 nearest neighbours of a point without one of its own                 RAD_ALL.25351 
     &     FILLIN,                                                                 RAD_ALL.25352 
C     !  Used to fill in points without a clearly-defined tropopause               RAD_ALL.25353 
     &     DTI              ! Default tropopause index for points where            RAD_ALL.25354 
C     ! not even one nearest neighbour has a well-defined tropopause               RAD_ALL.25355 
                                                                                   RAD_ALL.25356 
      REAL PJP1, PJ, PJM1,  ! Pressures at half levels J+1/J/J-1                   RAD_ALL.25357 
     &     P_EXNER_FULL_J,  ! Exner pressures at centres of current                RAD_ALL.25358 
     &     P_EXNER_FULL_JM1,!                      layer & one below               RAD_ALL.25359 
     &     DEL_EXNER_J,     ! Differences of Exner pressure across                 RAD_ALL.25360 
     &     DEL_EXNER_JM1,   !                           half layers                RAD_ALL.25361 
     &     DENOM            ! Denominator in lapse rate expression                 RAD_ALL.25362 
C*L------------------COMDECK C_G----------------------------------------           C_G.2     
C G IS MEAN ACCEL DUE TO GRAVITY AT EARTH'S SURFACE                                C_G.3     
      REAL G                                                                       C_G.4     
                                                                                   C_G.5     
      PARAMETER(G=9.80665)                                                         C_G.6     
C*----------------------------------------------------------------------           C_G.7     
C*L------------------COMDECK C_R_CP-------------------------------------           C_R_CP.2     
C R IS GAS CONSTANT FOR DRY AIR                                                    C_R_CP.3     
C CP IS SPECIFIC HEAT OF DRY AIR AT CONSTANT PRESSURE                              C_R_CP.4     
C PREF IS REFERENCE SURFACE PRESSURE                                               C_R_CP.5     
      REAL R,CP,KAPPA,PREF  ! PREFTOKI                                             C_R_CP.6     
                                                                                   C_R_CP.7     
      PARAMETER(R=287.05,                                                          C_R_CP.8     
     &          CP=1005.,                                                          C_R_CP.9     
     &          KAPPA=R/CP,                                                        C_R_CP.10    
     &          PREF=100000.)                                                      C_R_CP.11    
C*----------------------------------------------------------------------           C_R_CP.12    
                                                                                   C_R_CP.13    
C*L------------------COMDECK C_LAPSE ----------------------------------            C_LAPSE.2     
      REAL LAPSE,LAPSE_TROP                                                        C_LAPSE.3     
      PARAMETER(LAPSE=0.0065)     !  NEAR SURFACE LAPSE RATE                       C_LAPSE.4     
      PARAMETER(LAPSE_TROP=0.002) !  TROPOPAUSE LAPSE RATE                         C_LAPSE.5     
C*----------------------------------------------------------------------           C_LAPSE.6     
                                                                                   RAD_ALL.25366 
      REAL CP_OVER_G, P_EXNER_500, P_EXNER_50                                      RAD_ALL.25367 
      PARAMETER ( CP_OVER_G = CP / G )                                             RAD_ALL.25368 
                                                                                   RAD_ALL.25369 
C----------------------------------------------------------------------            RAD_ALL.25370 
C*L------------------ COMDECK P_EXNERC ---------------------------------           P_EXNERC.2     
C statement function to define exner pressure at model levels                      P_EXNERC.3     
C from exner pressures at model half-levels                                        P_EXNERC.4     
      REAL P_EXNER_C                                                               P_EXNERC.5     
      REAL R_P_EXNER_C                                                             P_EXNERC.6     
      REAL P_EXU_DUM,P_EXL_DUM,PU_DUM,PL_DUM,KAPPA_DUM                             P_EXNERC.7     
C consistent with geopotential see eqn 26 DOC PAPER 10                             P_EXNERC.15    
      P_EXNER_C(P_EXU_DUM,P_EXL_DUM,PU_DUM,PL_DUM,KAPPA_DUM) =                     P_EXNERC.16    
     & (P_EXU_DUM*PU_DUM - P_EXL_DUM*PL_DUM)/                                      P_EXNERC.17    
     & ( (PU_DUM-PL_DUM)*(KAPPA_DUM + 1) )                                         P_EXNERC.18    
      R_P_EXNER_C(P_EXU_DUM,P_EXL_DUM,PU_DUM,PL_DUM,KAPPA_DUM) =                   P_EXNERC.19    
     & ( (PU_DUM-PL_DUM)*(KAPPA_DUM + 1) )/                                        P_EXNERC.20    
     & (P_EXU_DUM*PU_DUM - P_EXL_DUM*PL_DUM)                                       P_EXNERC.21    
                                                                                   P_EXNERC.23    
C*------------------- --------------------------------------------------           P_EXNERC.24    
                                                                                   P_EXNERC.25    
                                                                                   RAD_ALL.25372 
C     ! 1.  Set up local constants and initialise arrays                           RAD_ALL.25373 
                                                                                   RAD_ALL.25374 
      P_EXNER_500 = (500./1000.)**KAPPA                                            RAD_ALL.25375 
      P_EXNER_50  =  (50./1000.)**KAPPA                                            RAD_ALL.25376 
                                                                                   RAD_ALL.25377 
      DTI = ( MIN_TROP_LEVEL + MAX_TROP_LEVEL ) / 2                                RAD_ALL.25378 
                                                                                   RAD_ALL.25379 
      DO I=1, L2                                                                   RAD_ALL.25380 
        LTROP(I) = .TRUE.                                                          RAD_ALL.25381 
      ENDDO                                                                        RAD_ALL.25382 
                                                                                   RAD_ALL.25383 
CL    ! Compute lapse rate between full levels: equation 3.16, UMDP S1             RAD_ALL.25384 
                                                                                   RAD_ALL.25385 
      DO J=MIN_TROP_LEVEL+1, MIN(MAX_TROP_LEVEL+1,P_LEVELS)                        RAD_ALL.25386 
        DO I=1, L2                                                                 RAD_ALL.25387 
C         ! Exner pressure at full levels                                          RAD_ALL.25388 
          PJP1 = AKH(J+1) + BKH(J+1) * PSTAR(I)                                    RAD_ALL.25389 
          PJ   = AKH(J)   + BKH(J)   * PSTAR(I)                                    RAD_ALL.25390 
          PJM1 = AKH(J-1) + BKH(J-1) * PSTAR(I)                                    RAD_ALL.25391 
          P_EXNER_FULL_J   = P_EXNER_C                                             RAD_ALL.25392 
     &      (P_EXNER_HALF(I,J+1), P_EXNER_HALF(I,J),  PJP1, PJ, KAPPA)             RAD_ALL.25393 
          P_EXNER_FULL_JM1 = P_EXNER_C                                             RAD_ALL.25394 
     &      (P_EXNER_HALF(I,J),  P_EXNER_HALF(I,J-1), PJ, PJM1, KAPPA)             RAD_ALL.25395 
C         ! Exner pressure difference across half layers                           RAD_ALL.25396 
          DEL_EXNER_J   = P_EXNER_HALF(I,J) - P_EXNER_FULL_J                       RAD_ALL.25397 
          DEL_EXNER_JM1 = P_EXNER_FULL_JM1  - P_EXNER_HALF(I,J)                    RAD_ALL.25398 
C         ! Denominator                                                            RAD_ALL.25399 
          DENOM = ( T(I,J-1) * DEL_EXNER_JM1 / P_EXNER_FULL_JM1                    RAD_ALL.25400 
     &                + T(I,J) * DEL_EXNER_J / P_EXNER_FULL_J )                    RAD_ALL.25401 
C         !  Lapse rate between level j-1 and j                                    RAD_ALL.25402 
          LAPSE_RATE(I,J) = ( T(I,J-1) - T(I,J) )/( CP_OVER_G * DENOM )            RAD_ALL.25403 
        ENDDO                                                                      RAD_ALL.25404 
      ENDDO                                                                        RAD_ALL.25405 
                                                                                   RAD_ALL.25406 
CL    ! 2.  Find level of tropopause, where it is well defined                     RAD_ALL.25407 
                                                                                   RAD_ALL.25408 
      DO J=MIN_TROP_LEVEL+1, MAX_TROP_LEVEL                                        RAD_ALL.25409 
                                                                                   RAD_ALL.25410 
C 'J+1' level for lapse rate test; allows J iteration up to P_LEVELS               RAD_ALL.25411 
        JP1=MIN(J+1,P_LEVELS)                                                      RAD_ALL.25412 
                                                                                   RAD_ALL.25413 
        DO I=1, L2                                                                 RAD_ALL.25414 
C         ! Exner pressure at full levels                                          RAD_ALL.25415 
          PJP1 = AKH(J+1) + BKH(J+1) * PSTAR(I)                                    RAD_ALL.25416 
          PJ   = AKH(J)   + BKH(J)   * PSTAR(I)                                    RAD_ALL.25417 
          PJM1 = AKH(J-1) + BKH(J-1) * PSTAR(I)                                    RAD_ALL.25418 
          P_EXNER_FULL_J = P_EXNER_C                                               RAD_ALL.25419 
     &       (P_EXNER_HALF(I,J+1), P_EXNER_HALF(I,J), PJP1, PJ, KAPPA)             RAD_ALL.25420 
          P_EXNER_FULL_JM1 = P_EXNER_C                                             RAD_ALL.25421 
     &       (P_EXNER_HALF(I,J), P_EXNER_HALF(I,J-1), PJ, PJM1, KAPPA)             RAD_ALL.25422 
                                                                                   RAD_ALL.25423 
C         ! Not-quite-WMO criteria for interval containing tropopause              RAD_ALL.25424 
C         ! (where 'interval' stretches between layer centres j and j-1)           RAD_ALL.25425 
                                                                                   RAD_ALL.25426 
          IF (  P_EXNER_FULL_JM1 .GT. P_EXNER_50  .AND.                            RAD_ALL.25427 
     &            P_EXNER_FULL_J .LT. P_EXNER_500 .AND.                            RAD_ALL.25428 
     &           LAPSE_RATE(I,J) .LT. LAPSE_TROP  .AND.                            RAD_ALL.25429 
     &         LAPSE_RATE(I,JP1) .LT. LAPSE_TROP  .AND. LTROP(I) )                 RAD_ALL.25430 
     &    THEN                                                                     RAD_ALL.25431 
            LTROP(I)=.FALSE.                                                       RAD_ALL.25432 
            IT(I) = J                                                              RAD_ALL.25433 
          ENDIF                                                                    RAD_ALL.25434 
        ENDDO                                                                      RAD_ALL.25435 
      ENDDO                                                                        RAD_ALL.25436 
                                                                                   RAD_ALL.25437 
                                                                                   RAD_ALL.25438 
CL    ! 4.  Fill in where the above criteria did not find a tropopause             RAD_ALL.25439 
                                                                                   RAD_ALL.25440 
C     !  Run through all internal points and where no tropopause was               RAD_ALL.25441 
C     !   found, set the level to the average of those found at the 8              RAD_ALL.25442 
C     !   surrounding points.  If none of these did find one, then DTI,            RAD_ALL.25443 
C     !   the middle of the permitted range, is set.                               RAD_ALL.25444 
C     !   (Using integers means rounding down - probably the best thing            RAD_ALL.25445 
C     !   overall as the lack of vertical resolution to pick out a                 RAD_ALL.25446 
C     !   tropopause precisely is likely to mean they'll be diagnosed              RAD_ALL.25447 
C     !   too high, if anything, at neighbouring points.  Similarly                RAD_ALL.25448 
C     !   it's not worth worrying about level number being non-linear              RAD_ALL.25449 
C     !   in height or pressure at those points where ipso facto the               RAD_ALL.25450 
C     !   result is so arbitrary.)                                                 RAD_ALL.25451 
                                                                                   RAD_ALL.25452 
      DO J=1, L2/ROW_LENGTH-2                                                      RAD_ALL.25453 
        DO I=2, ROW_LENGTH-1                                                       RAD_ALL.25454 
          POINT = I+J*ROW_LENGTH                                                   RAD_ALL.25455 
          IF ( LTROP(POINT) ) THEN                                                 RAD_ALL.25456 
            FILLIN = 0                                                             RAD_ALL.25457 
            NNEIGH = 0                                                             RAD_ALL.25458 
            IF ( .NOT. LTROP(POINT-1-ROW_LENGTH) ) THEN                            RAD_ALL.25459 
              FILLIN = FILLIN + IT(POINT-1-ROW_LENGTH)                             RAD_ALL.25460 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25461 
            ENDIF                                                                  RAD_ALL.25462 
            IF ( .NOT. LTROP(POINT+1-ROW_LENGTH) ) THEN                            RAD_ALL.25463 
              FILLIN = FILLIN + IT(POINT+1-ROW_LENGTH)                             RAD_ALL.25464 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25465 
            ENDIF                                                                  RAD_ALL.25466 
            IF ( .NOT. LTROP(POINT-1+ROW_LENGTH) ) THEN                            RAD_ALL.25467 
              FILLIN = FILLIN + IT(POINT-1+ROW_LENGTH)                             RAD_ALL.25468 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25469 
            ENDIF                                                                  RAD_ALL.25470 
            IF ( .NOT. LTROP(POINT+1+ROW_LENGTH) ) THEN                            RAD_ALL.25471 
              FILLIN = FILLIN + IT(POINT+1+ROW_LENGTH)                             RAD_ALL.25472 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25473 
            ENDIF                                                                  RAD_ALL.25474 
            IF ( .NOT. LTROP(POINT-1) ) THEN                                       RAD_ALL.25475 
              FILLIN = FILLIN + IT(POINT-1)                                        RAD_ALL.25476 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25477 
            ENDIF                                                                  RAD_ALL.25478 
            IF ( .NOT. LTROP(POINT+1) ) THEN                                       RAD_ALL.25479 
              FILLIN = FILLIN + IT(POINT+1)                                        RAD_ALL.25480 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25481 
            ENDIF                                                                  RAD_ALL.25482 
            IF ( .NOT. LTROP(POINT-ROW_LENGTH) ) THEN                              RAD_ALL.25483 
              FILLIN = FILLIN + IT(POINT-ROW_LENGTH)                               RAD_ALL.25484 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25485 
            ENDIF                                                                  RAD_ALL.25486 
            IF ( .NOT. LTROP(POINT+ROW_LENGTH) ) THEN                              RAD_ALL.25487 
              FILLIN = FILLIN + IT(POINT+ROW_LENGTH)                               RAD_ALL.25488 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25489 
            ENDIF                                                                  RAD_ALL.25490 
            IF ( NNEIGH .EQ. 0 ) THEN                                              RAD_ALL.25491 
               IT(POINT) = DTI                                                     RAD_ALL.25492 
             ELSE                                                                  RAD_ALL.25493 
               IT(POINT) = FILLIN / NNEIGH                                         RAD_ALL.25494 
            ENDIF                                                                  RAD_ALL.25495 
          ENDIF                                                                    RAD_ALL.25496 
        ENDDO                                                                      RAD_ALL.25497 
      ENDDO                                                                        RAD_ALL.25498 
                                                                                   RAD_ALL.25499 
C     !  If the grid wraps round, the first & last points on each row              RAD_ALL.25500 
C     !    are in fact internal points & can be filled in similarly if             RAD_ALL.25501 
C     !    a few indices are altered.  (Indeed, they must be if the                RAD_ALL.25502 
C     !    model is to give the same answers regardless of                         RAD_ALL.25503 
C     !    decomposition - there may be no physical justification for              RAD_ALL.25504 
C     !    being so scrupulous.)                                                   RAD_ALL.25505 
C     !  If not, again such points return DTI if no tropopause was                 RAD_ALL.25506 
C     !    found on the WMO criteria.                                              RAD_ALL.25507 
                                                                                   RAD_ALL.25508 
      IF ( WRAP ) THEN                                                             RAD_ALL.25509 
         DO J=1, L2/ROW_LENGTH-2                                                   RAD_ALL.25510 
           POINT = 1+J*ROW_LENGTH                                                  RAD_ALL.25511 
           IF ( LTROP(POINT) ) THEN                                                RAD_ALL.25512 
             FILLIN = 0                                                            RAD_ALL.25513 
             NNEIGH = 0                                                            RAD_ALL.25514 
             IF ( .NOT. LTROP(POINT-1) ) THEN                                      RAD_ALL.25515 
               FILLIN = FILLIN + IT(POINT-1)                                       RAD_ALL.25516 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25517 
             ENDIF                                                                 RAD_ALL.25518 
             IF ( .NOT. LTROP(POINT+1-ROW_LENGTH) ) THEN                           RAD_ALL.25519 
               FILLIN = FILLIN + IT(POINT+1-ROW_LENGTH)                            RAD_ALL.25520 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25521 
             ENDIF                                                                 RAD_ALL.25522 
             IF ( .NOT. LTROP(POINT-1+2*ROW_LENGTH) ) THEN                         RAD_ALL.25523 
               FILLIN = FILLIN + IT(POINT-1+2*ROW_LENGTH)                          RAD_ALL.25524 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25525 
             ENDIF                                                                 RAD_ALL.25526 
             IF ( .NOT. LTROP(POINT+1+ROW_LENGTH) ) THEN                           RAD_ALL.25527 
               FILLIN = FILLIN + IT(POINT+1+ROW_LENGTH)                            RAD_ALL.25528 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25529 
             ENDIF                                                                 RAD_ALL.25530 
             IF ( .NOT. LTROP(POINT-1+ROW_LENGTH) ) THEN                           RAD_ALL.25531 
               FILLIN = FILLIN + IT(POINT-1+ROW_LENGTH)                            RAD_ALL.25532 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25533 
             ENDIF                                                                 RAD_ALL.25534 
             IF ( .NOT. LTROP(POINT+1) ) THEN                                      RAD_ALL.25535 
               FILLIN = FILLIN + IT(POINT+1)                                       RAD_ALL.25536 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25537 
             ENDIF                                                                 RAD_ALL.25538 
             IF ( .NOT. LTROP(POINT-ROW_LENGTH) ) THEN                             RAD_ALL.25539 
               FILLIN = FILLIN + IT(POINT-ROW_LENGTH)                              RAD_ALL.25540 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25541 
             ENDIF                                                                 RAD_ALL.25542 
             IF ( .NOT. LTROP(POINT+ROW_LENGTH) ) THEN                             RAD_ALL.25543 
               FILLIN = FILLIN + IT(POINT+ROW_LENGTH)                              RAD_ALL.25544 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25545 
             ENDIF                                                                 RAD_ALL.25546 
             IF ( NNEIGH .EQ. 0 ) THEN                                             RAD_ALL.25547 
                IT(POINT) = DTI                                                    RAD_ALL.25548 
              ELSE                                                                 RAD_ALL.25549 
                IT(POINT) = FILLIN / NNEIGH                                        RAD_ALL.25550 
             ENDIF                                                                 RAD_ALL.25551 
           ENDIF                                                                   RAD_ALL.25552 
           POINT = (1+J)*ROW_LENGTH                                                RAD_ALL.25553 
           IF ( LTROP(POINT) ) THEN                                                RAD_ALL.25554 
             FILLIN = 0                                                            RAD_ALL.25555 
             NNEIGH = 0                                                            RAD_ALL.25556 
             IF ( .NOT. LTROP(POINT-1-ROW_LENGTH) ) THEN                           RAD_ALL.25557 
               FILLIN = FILLIN + IT(POINT-1-ROW_LENGTH)                            RAD_ALL.25558 
              NNEIGH = NNEIGH + 1                                                  RAD_ALL.25559 
             ENDIF                                                                 RAD_ALL.25560 
             IF ( .NOT. LTROP(POINT+1-2*ROW_LENGTH) ) THEN                         RAD_ALL.25561 
               FILLIN = FILLIN + IT(POINT+1-2*ROW_LENGTH)                          RAD_ALL.25562 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25563 
             ENDIF                                                                 RAD_ALL.25564 
             IF ( .NOT. LTROP(POINT-1+ROW_LENGTH) ) THEN                           RAD_ALL.25565 
               FILLIN = FILLIN + IT(POINT-1+ROW_LENGTH)                            RAD_ALL.25566 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25567 
             ENDIF                                                                 RAD_ALL.25568 
             IF ( .NOT. LTROP(POINT+1) ) THEN                                      RAD_ALL.25569 
               FILLIN = FILLIN + IT(POINT+1)                                       RAD_ALL.25570 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25571 
             ENDIF                                                                 RAD_ALL.25572 
             IF ( .NOT. LTROP(POINT-1) ) THEN                                      RAD_ALL.25573 
               FILLIN = FILLIN + IT(POINT-1)                                       RAD_ALL.25574 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25575 
             ENDIF                                                                 RAD_ALL.25576 
             IF ( .NOT. LTROP(POINT+1-ROW_LENGTH) ) THEN                           RAD_ALL.25577 
               FILLIN = FILLIN + IT(POINT+1-ROW_LENGTH)                            RAD_ALL.25578 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25579 
             ENDIF                                                                 RAD_ALL.25580 
             IF ( .NOT. LTROP(POINT-ROW_LENGTH) ) THEN                             RAD_ALL.25581 
               FILLIN = FILLIN + IT(POINT-ROW_LENGTH)                              RAD_ALL.25582 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25583 
             ENDIF                                                                 RAD_ALL.25584 
             IF ( .NOT. LTROP(POINT+ROW_LENGTH) ) THEN                             RAD_ALL.25585 
               FILLIN = FILLIN + IT(POINT+ROW_LENGTH)                              RAD_ALL.25586 
               NNEIGH = NNEIGH + 1                                                 RAD_ALL.25587 
             ENDIF                                                                 RAD_ALL.25588 
             IF ( NNEIGH .EQ. 0 ) THEN                                             RAD_ALL.25589 
                IT(POINT) = DTI                                                    RAD_ALL.25590 
              ELSE                                                                 RAD_ALL.25591 
                IT(POINT) = FILLIN / NNEIGH                                        RAD_ALL.25592 
             ENDIF                                                                 RAD_ALL.25593 
           ENDIF                                                                   RAD_ALL.25594 
         ENDDO                                                                     RAD_ALL.25595 
       ELSE ! if not WRAP                                                          RAD_ALL.25596 
         DO J=1, L2/ROW_LENGTH-2                                                   RAD_ALL.25597 
           IF ( LTROP(1+J*ROW_LENGTH) )   IT(1+J*ROW_LENGTH)   = DTI               RAD_ALL.25598 
           IF ( LTROP((1+J)*ROW_LENGTH) ) IT((1+J)*ROW_LENGTH) = DTI               RAD_ALL.25599 
         ENDDO                                                                     RAD_ALL.25600 
      ENDIF                                                                        RAD_ALL.25601 
C                                                                                  RAD_ALL.25602 
      DO J=1, ROW_LENGTH                                                           RAD_ALL.25603 
        IF ( LTROP(J) )      IT(J)      = DTI                                      RAD_ALL.25604 
        IF ( LTROP(L2+1-J) ) IT(L2+1-J) = DTI                                      RAD_ALL.25605 
      ENDDO                                                                        RAD_ALL.25606 
C                                                                                  RAD_ALL.25607 
      RETURN                                                                       RAD_ALL.25608 
      END 
