Back to home page

MITgcm

 
 

    


File indexing completed on 2023-07-14 05:10:23 UTC

view on githubraw file Latest commit de57a2ec on 2023-07-13 16:55:13 UTC
8f7d13d0c9 Jean*0001 #include "ECCO_OPTIONS.h"
f1e8c1d7ee Gael*0002 
                0003       subroutine cost_gencost_sstv4(
f8e779c983 antn*0004      I                     myThid
f1e8c1d7ee Gael*0005      &                   )
                0006 
                0007 c     ==================================================================
                0008 c     SUBROUTINE cost_gencost_sstv4
                0009 c     ==================================================================
                0010 c
                0011 c     o Evaluate cost function contributions of sea surface temperature.
                0012 c       (Daily Pointwise and then Large Scale)
                0013 c
                0014 c       started: Gael Forget, Oct-2009
                0015 c
                0016 c     ==================================================================
                0017 c     SUBROUTINE cost_gencost_sstv4
                0018 c     ==================================================================
                0019 
                0020       implicit none
                0021 
                0022 c     == global variables ==
                0023 
                0024 #include "EEPARAMS.h"
                0025 #include "SIZE.h"
                0026 #include "PARAMS.h"
                0027 #include "GRID.h"
f09238ab8f Gael*0028 #ifdef ALLOW_CAL
                0029 # include "cal.h"
                0030 #endif
                0031 #ifdef ALLOW_ECCO
13d362b8c1 Ou W*0032 # include "ECCO_SIZE.h"
                0033 # include "ECCO.h"
f09238ab8f Gael*0034 #endif
ebac42f582 An T*0035 #ifdef ALLOW_SMOOTH
f09238ab8f Gael*0036 # include "SMOOTH.h"
ebac42f582 An T*0037 #endif
f1e8c1d7ee Gael*0038 
                0039 c     == routine arguments ==
                0040 
f8e779c983 antn*0041       integer myThid
f1e8c1d7ee Gael*0042 
                0043 #ifdef ALLOW_SMOOTH
                0044 #ifdef ALLOW_GENCOST_CONTRIBUTION
                0045 c     == local variables ==
                0046 
                0047       integer bi,bj
                0048       integer i,j,k
                0049       integer itlo,ithi
                0050       integer jtlo,jthi
                0051       integer jmin,jmax
                0052       integer imin,imax
                0053       integer irec,jrec,krec
                0054       integer ilps
                0055 
                0056       logical doglobalread
                0057       logical ladinit
                0058 
f09238ab8f Gael*0059       _RL mydummy
f8e779c983 antn*0060       _RL mybar(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0061       _RL anom_sst(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0062       _RL obs_sst (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0063       _RL nb_sst (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0064       _RL msk_sst (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0065       _RL tmp_sst (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
f1e8c1d7ee Gael*0066       _RL spval
                0067 
                0068       _RL junk,junkweight
                0069 
                0070       integer ndaysave
                0071       _RL ndaysaveRL
                0072 
3248fb8b5a Gael*0073       integer k2, k2_lsc
                0074 
de57a2ec4b Mart*0075       character*(MAX_LEN_FNAM) fname
                0076       character*(MAX_LEN_FNAM) fname2
e5310f7c13 Gael*0077 
                0078 #ifdef ALLOW_ECCO_DEBUG
                0079       character*(MAX_LEN_MBUF) msgBuf
                0080       INTEGER ioUnit
                0081 #endif
                0082       logical exst
f1e8c1d7ee Gael*0083 
                0084       _RL daytime
                0085       _RL diffsecs
                0086       integer il, localrec
                0087       integer dayiter
                0088       integer daydate(4)
                0089       integer difftime(4)
e36a586433 Jean*0090       integer tempDate_1
f1e8c1d7ee Gael*0091       integer middate(4)
e5310f7c13 Gael*0092       integer locstartdate(4)
f1e8c1d7ee Gael*0093       integer yday, ymod
                0094       integer md, dd, sd, ld, wd
                0095 
f09238ab8f Gael*0096       integer kgen, kgen_lsc
b9b1c4d04c Timo*0097       logical dosumsq
f1e8c1d7ee Gael*0098 
                0099 c     == external functions ==
                0100 
                0101       integer  ilnblnk
                0102       external ilnblnk
                0103 
                0104 c     == end of interface ==
                0105 
f8e779c983 antn*0106       jtlo = myByLo(myThid)
                0107       jthi = myByHi(myThid)
                0108       itlo = myBxLo(myThid)
                0109       ithi = myBxHi(myThid)
f1e8c1d7ee Gael*0110       jmin = 1
f8e779c983 antn*0111       jmax = sNy
f1e8c1d7ee Gael*0112       imin = 1
f8e779c983 antn*0113       imax = sNx
f1e8c1d7ee Gael*0114 
e5310f7c13 Gael*0115 #ifdef ALLOW_ECCO_DEBUG
                0116       ioUnit=standardMessageUnit
                0117 #endif
                0118 
f1e8c1d7ee Gael*0119 c-- detect the relevant gencost indices
f09238ab8f Gael*0120       kgen=0
                0121       kgen_lsc=0
3248fb8b5a Gael*0122       k2_lsc=0
f1e8c1d7ee Gael*0123       do k=1,NGENCOST
3248fb8b5a Gael*0124        if (gencost_name(k).EQ.'sstv4-amsre') kgen=k
                0125        if (gencost_name(k).EQ.'sstv4-amsre-lsc') then
                0126         kgen_lsc=k
                0127         do k2 = 1, NGENPPROC
                0128          if (gencost_posproc(k2,kgen_lsc).EQ.'smooth') k2_lsc=k2
                0129         enddo
                0130        endif
f1e8c1d7ee Gael*0131       enddo
                0132 
f09238ab8f Gael*0133       if (kgen.NE.0) then
f8e779c983 antn*0134 c ------
f09238ab8f Gael*0135 
b9b1c4d04c Timo*0136       dosumsq=.TRUE.
11c3150c71 Mart*0137       call ecco_zero( gencost_weight(1-OLx,1-OLy,1,1,kgen),
                0138      &                1, zeroRL, myThid )
                0139       call ecco_zero( gencost_weight(1-OLx,1-OLy,1,1,kgen_lsc),
                0140      &                1, zeroRL, myThid )
cbd85e4123 Gael*0141       if ( gencost_errfile(kgen) .NE. ' ' )
11c3150c71 Mart*0142      &   call ecco_readwei( gencost_errfile(kgen),
                0143      &     gencost_weight(1-OLx,1-OLy,1,1,kgen),
                0144      &     1, 1, 1, dosumsq, myThid )
cbd85e4123 Gael*0145       if ( gencost_errfile(kgen_lsc) .NE. ' ' )
11c3150c71 Mart*0146      &   call ecco_readwei( gencost_errfile(kgen_lsc),
                0147      &     gencost_weight(1-OLx,1-OLy,1,1,kgen_lsc),
                0148      &     1, 1, 1, dosumsq, myThid )
cbd85e4123 Gael*0149 
f8e779c983 antn*0150       call cal_FullDate(19920101,0,locstartdate,myThid)
f09238ab8f Gael*0151 
f1e8c1d7ee Gael*0152 c--   First, read tiled data.
                0153       doglobalread = .false.
                0154       ladinit      = .false.
                0155 
f09238ab8f Gael*0156       ilps=ilnblnk( gencost_barfile(kgen) )
de57a2ec4b Mart*0157       write(fname,'(2a,i10.10)')
f09238ab8f Gael*0158      &     gencost_barfile(kgen)(1:ilps),'.',eccoiter
f1e8c1d7ee Gael*0159 
f09238ab8f Gael*0160       spval = gencost_spmin(kgen)
                0161       mydummy = gencost_dummy(kgen)
f8e779c983 antn*0162 
f1e8c1d7ee Gael*0163 cgf =======================================================
                0164 cgf PART 1: compute smooth SST cost term
                0165 cgf =======================================================
                0166 
                0167       ndaysave=7
                0168       ndaysaveRL=ndaysave
                0169 
                0170       do irec = 1, ndaysrec-ndaysave+1, 7
                0171 
                0172          do bj = jtlo,jthi
                0173           do bi = itlo,ithi
                0174            do j = jmin,jmax
                0175             do i = imin,imax
                0176               anom_sst(i,j,bi,bj)  = 0. _d 0
                0177               obs_sst(i,j,bi,bj)  = 0. _d 0
                0178               nb_sst(i,j,bi,bj)  = 0. _d 0
                0179               msk_sst(i,j,bi,bj)  = 0. _d 0
                0180             enddo
                0181            enddo
                0182           enddo
                0183          enddo
                0184 
                0185 c PART 1.1: compute running sample average over ndaysave
                0186 c ------------------------------------------------------
                0187 
                0188          do jrec=1,ndaysave
                0189 
                0190            krec=irec+jrec-1
                0191 
                0192 c get modeled sst:
101f75e5cd Gael*0193 #ifdef ALLOW_AUTODIFF
f09238ab8f Gael*0194              call active_read_xy( fname, mybar, krec, doglobalread,
f8e779c983 antn*0195      &                       ladinit, eccoiter, myThid,
f09238ab8f Gael*0196      &                       mydummy )
101f75e5cd Gael*0197 #else
                0198              CALL READ_REC_XY_RL( fname, mybar, kRec, 1, myThid )
                0199 #endif
f1e8c1d7ee Gael*0200 
                0201 c get observed sst:
e5310f7c13 Gael*0202              daytime = FLOAT(secondsperday*(krec-1)) + modelstart
                0203              dayiter = hoursperday*(krec-1)+modeliter0
f8e779c983 antn*0204              call cal_getdate( dayiter, daytime, daydate, myThid )
                0205              call cal_convdate( daydate,yday,md,dd,sd,ld,wd,myThid )
e5310f7c13 Gael*0206              ymod = locstartdate(1)/10000
                0207              if ( ymod .GE. yday ) then
e244d7989d Gael*0208                middate(1)=1
f8e779c983 antn*0209                call cal_FullDate(locstartdate(1),0,middate,myThid)
f1e8c1d7ee Gael*0210              else
e244d7989d Gael*0211                middate(1)=1
e36a586433 Jean*0212                tempDate_1 = yday*10000+100+1
f8e779c983 antn*0213                call cal_FullDate( tempDate_1, 0, middate, myThid)
f1e8c1d7ee Gael*0214              endif
f8e779c983 antn*0215              call cal_TimePassed( middate, daydate, difftime, myThid )
                0216              call cal_ToSeconds( difftime, diffsecs, myThid )
e5310f7c13 Gael*0217 c             localrec = floor(diffsecs/86400.) + 1
f1e8c1d7ee Gael*0218              localrec = int(diffsecs/86400.) + 1
                0219 
f09238ab8f Gael*0220              il=ilnblnk(gencost_datafile(kgen))
de57a2ec4b Mart*0221              write(fname2,'(2a,i4)')
f09238ab8f Gael*0222      &         gencost_datafile(kgen)(1:il), '_', yday
e5310f7c13 Gael*0223              inquire( file=fname2, exist=exst )
                0224 
                0225 #ifdef ALLOW_ECCO_DEBUG
                0226         WRITE(msgBuf,'(A,I4,A,I4,A,I10,A,1PE15.2)') 'sstv4 reading ',
                0227      &      yday,' ',ymod,' ',localrec,' ',diffsecs
                0228         CALL PRINT_MESSAGE( msgBuf, ioUnit, SQUEEZE_RIGHT, myThid )
                0229 c
f8e779c983 antn*0230         CALL CAL_PRINTDATE(middate,myThid)
                0231         CALL CAL_PRINTDATE(daydate,myThid)
                0232         CALL CAL_PRINTDATE(difftime,myThid)
e36a586433 Jean*0233 #endif
f1e8c1d7ee Gael*0234 
e244d7989d Gael*0235              if ( ( localrec .GT. 0 ).AND.(diffsecs .GT. 0.d0) ) then
f8e779c983 antn*0236                CALL READ_REC_3D_RL( fname2, cost_iprec, 1,
                0237      &                              tmp_sst, localrec, 1, myThid )
f1e8c1d7ee Gael*0238              else
                0239               do bj = jtlo,jthi
                0240                 do bi = itlo,ithi
                0241                     do j = jmin,jmax
                0242                       do i = imin,imax
                0243                          tmp_sst(i,j,bi,bj) = spval
                0244                       enddo
                0245                     enddo
                0246                 enddo
                0247               enddo
                0248              endif
                0249 
                0250 c accumulate obs and misfit:
                0251              do bj = jtlo,jthi
                0252               do bi = itlo,ithi
                0253                do j = jmin,jmax
                0254                 do i = imin,imax
                0255                  if ( (tmp_sst(i,j,bi,bj).GT.spval).AND.
                0256      &                (maskc(i,j,1,bi,bj).EQ.1.) ) then
                0257                   anom_sst(i,j,bi,bj)= anom_sst(i,j,bi,bj)+
f09238ab8f Gael*0258      &               mybar(i,j,bi,bj)-tmp_sst(i,j,bi,bj)
f1e8c1d7ee Gael*0259                   obs_sst(i,j,bi,bj)= obs_sst(i,j,bi,bj)+
e36a586433 Jean*0260      &               tmp_sst(i,j,bi,bj)
f1e8c1d7ee Gael*0261                   nb_sst(i,j,bi,bj)=nb_sst(i,j,bi,bj)+1. _d 0
                0262                  endif
                0263                 enddo
                0264                enddo
                0265               enddo
                0266              enddo
                0267 
                0268          enddo !do jrec=1,ndaysave
                0269 
                0270 c average obs and misfit:
                0271          do bj = jtlo,jthi
                0272           do bi = itlo,ithi
                0273            do j = jmin,jmax
                0274             do i = imin,imax
                0275              if ( nb_sst(i,j,bi,bj) .NE. 0. ) then
e36a586433 Jean*0276               obs_sst(i,j,bi,bj) =
f1e8c1d7ee Gael*0277      &            obs_sst(i,j,bi,bj)/nb_sst(i,j,bi,bj)
e36a586433 Jean*0278               anom_sst(i,j,bi,bj) =
f1e8c1d7ee Gael*0279      &            anom_sst(i,j,bi,bj)/nb_sst(i,j,bi,bj)
                0280               msk_sst(i,j,bi,bj) = 1. _d 0
                0281              endif
                0282             enddo
                0283            enddo
                0284           enddo
                0285          enddo
                0286 
                0287 c PART 1.2: smooth anom_sst in space
                0288 c ----------------------------------------
                0289 
                0290 #ifdef ALLOW_GENCOST_SSTV4_OUTPUT
de57a2ec4b Mart*0291          write(fname2,'(1a)') 'sstdiff_raw'
f8e779c983 antn*0292          CALL WRITE_REC_3D_RL( fname2, precFloat32, 1,
                0293      &                         anom_sst, irec, 1, myThid )
de57a2ec4b Mart*0294          write(fname2,'(1a)') 'sstobs_raw'
f8e779c983 antn*0295          CALL WRITE_REC_3D_RL( fname2, precFloat32, 1,
                0296      &                         obs_sst, irec, 1, myThid )
f1e8c1d7ee Gael*0297 #endif
                0298 
3248fb8b5a Gael*0299          if ( useSMOOTH.AND.(k2_lsc.GT.0) )
7b8b86ab99 Timo*0300      &     call smooth_hetero2d(anom_sst,maskInC,
3248fb8b5a Gael*0301      &     gencost_posproc_c(k2_lsc,kgen_lsc),
f8e779c983 antn*0302      &     gencost_posproc_i(k2_lsc,kgen_lsc),myThid)
f1e8c1d7ee Gael*0303 
                0304 #ifdef ALLOW_GENCOST_SSTV4_OUTPUT
3248fb8b5a Gael*0305          if ( useSMOOTH.AND.(k2_lsc.GT.0) )
7b8b86ab99 Timo*0306      &     call smooth_hetero2d(obs_sst,maskInC,
3248fb8b5a Gael*0307      &     gencost_posproc_c(k2_lsc,kgen_lsc),
f8e779c983 antn*0308      &     gencost_posproc_i(k2_lsc,kgen_lsc),myThid)
f1e8c1d7ee Gael*0309 
de57a2ec4b Mart*0310          write(fname2,'(1a)') 'sstdiff_smooth'
f8e779c983 antn*0311          CALL WRITE_REC_3D_RL( fname2, precFloat32, 1,
                0312      &                         anom_sst, irec, 1, myThid )
de57a2ec4b Mart*0313          write(fname2,'(1a)') 'sstobs_smooth'
f8e779c983 antn*0314          CALL WRITE_REC_3D_RL( fname2, precFloat32, 1,
                0315      &                         obs_sst, irec, 1, myThid )
e36a586433 Jean*0316 #endif
f1e8c1d7ee Gael*0317 
                0318 c PART 1.3: compute cost function term
                0319 c ------------------------------------
                0320 
                0321          do bj = jtlo,jthi
                0322           do bi = itlo,ithi
                0323            do j = jmin,jmax
                0324             do i = imin,imax
                0325              junk = anom_sst(i,j,bi,bj)
f09238ab8f Gael*0326              junkweight = gencost_weight(i,j,bi,bj,kgen_lsc)*
f1e8c1d7ee Gael*0327      &          maskc(i,j,1,bi,bj)
f4cfad5b95 Ou W*0328              objf_gencost(bi,bj,kgen_lsc) =
                0329      &          objf_gencost(bi,bj,kgen_lsc)
f1e8c1d7ee Gael*0330      &            +junk*junk*junkweight/ndaysaveRL
                0331              if ( (junkweight.GT.0.).AND.(nb_sst(i,j,bi,bj).GT.0.) )
f4cfad5b95 Ou W*0332      &          num_gencost(bi,bj,kgen_lsc) =
                0333      &          num_gencost(bi,bj,kgen_lsc) + 1. _d 0 /ndaysaveRL
f1e8c1d7ee Gael*0334             enddo
                0335            enddo
                0336           enddo
                0337          enddo
                0338 
                0339       enddo
                0340 
                0341 cgf =======================================================
                0342 cgf PART 2: compute raw SST cost term
                0343 cgf =======================================================
                0344 
                0345       do irec = 1, ndaysrec
                0346 
                0347 c get modeled sst:
101f75e5cd Gael*0348 #ifdef ALLOW_AUTODIFF
f09238ab8f Gael*0349         call active_read_xy( fname, mybar, irec, doglobalread,
f8e779c983 antn*0350      &                       ladinit, eccoiter, myThid,
f09238ab8f Gael*0351      &                       mydummy )
101f75e5cd Gael*0352 #else
                0353         CALL READ_REC_XY_RL( fname, mybar, iRec, 1, myThid )
                0354 #endif
f1e8c1d7ee Gael*0355 
                0356 c get observed sst:
e5310f7c13 Gael*0357              daytime = FLOAT(secondsperday*(irec-1)) + modelstart
                0358              dayiter = hoursperday*(irec-1)+modeliter0
f8e779c983 antn*0359              call cal_getdate( dayiter, daytime, daydate, myThid )
                0360              call cal_convdate( daydate,yday,md,dd,sd,ld,wd,myThid )
e5310f7c13 Gael*0361              ymod = locstartdate(1)/10000
                0362              if ( ymod .GE. yday ) then
e244d7989d Gael*0363                middate(1)=1
f8e779c983 antn*0364                call cal_FullDate(locstartdate(1),0,middate,myThid)
f1e8c1d7ee Gael*0365              else
e244d7989d Gael*0366                middate(1)=1
                0367                tempDate_1 = yday*10000+100+1
f8e779c983 antn*0368                call cal_FullDate( tempDate_1, 0, middate, myThid)
f1e8c1d7ee Gael*0369              endif
f8e779c983 antn*0370              call cal_TimePassed( middate, daydate, difftime, myThid )
                0371              call cal_ToSeconds( difftime, diffsecs, myThid )
e5310f7c13 Gael*0372 c             localrec = floor(diffsecs/86400.) + 1
f1e8c1d7ee Gael*0373              localrec = int(diffsecs/86400.) + 1
                0374 
f09238ab8f Gael*0375              il=ilnblnk(gencost_datafile(kgen))
de57a2ec4b Mart*0376              write(fname2,'(2a,i4)')
f09238ab8f Gael*0377      &         gencost_datafile(kgen)(1:il), '_', yday
e5310f7c13 Gael*0378              inquire( file=fname2, exist=exst )
                0379 
                0380 #ifdef ALLOW_ECCO_DEBUG
                0381         WRITE(msgBuf,'(A,I4,A,I4,A,I10,A,1PE15.2)') 'sstv4 reading ',
                0382      &      yday,' ',ymod,' ',localrec,' ',diffsecs
                0383         CALL PRINT_MESSAGE( msgBuf, ioUnit, SQUEEZE_RIGHT, myThid )
                0384 c
f8e779c983 antn*0385         CALL CAL_PRINTDATE(middate,myThid)
                0386         CALL CAL_PRINTDATE(daydate,myThid)
                0387         CALL CAL_PRINTDATE(difftime,myThid)
e5310f7c13 Gael*0388 #endif
f1e8c1d7ee Gael*0389 
e244d7989d Gael*0390              if ( ( localrec .GT. 0 ).AND.(diffsecs .GT. 0.d0) ) then
f8e779c983 antn*0391                call READ_REC_3D_RL( fname2, cost_iprec, 1,
                0392      &                              tmp_sst, localrec, 1, myThid )
f1e8c1d7ee Gael*0393              else
                0394               do bj = jtlo,jthi
                0395                 do bi = itlo,ithi
                0396                     do j = jmin,jmax
                0397                       do i = imin,imax
                0398                          tmp_sst(i,j,bi,bj) = spval
                0399                       enddo
                0400                     enddo
                0401                 enddo
                0402               enddo
                0403              endif
                0404 
                0405 c compute misfit:
                0406          do bj = jtlo,jthi
                0407           do bi = itlo,ithi
                0408            do j = jmin,jmax
                0409             do i = imin,imax
                0410              if ( (tmp_sst(i,j,bi,bj).GT.spval).AND.
                0411      &            (maskc(i,j,1,bi,bj).EQ.1.) ) then
                0412               anom_sst(i,j,bi,bj) =
f09238ab8f Gael*0413      &               mybar(i,j,bi,bj)-tmp_sst(i,j,bi,bj)
f1e8c1d7ee Gael*0414               msk_sst(i,j,bi,bj) = 1. _d 0
                0415              else
                0416               anom_sst(i,j,bi,bj) = 0. _d 0
                0417               msk_sst(i,j,bi,bj) = 0. _d 0
                0418              endif
                0419             enddo
                0420            enddo
                0421           enddo
                0422          enddo
                0423 
                0424 #ifdef ALLOW_GENCOST_SSTV4_OUTPUT
de57a2ec4b Mart*0425          write(fname2,'(1a)') 'sstdiff_point'
f8e779c983 antn*0426          CALL WRITE_REC_3D_RL( fname2, precFloat32, 1,
                0427      &                         anom_sst, irec, 1, myThid )
e36a586433 Jean*0428 #endif
f1e8c1d7ee Gael*0429 
                0430 c compute cost:
                0431 
                0432          do bj = jtlo,jthi
                0433           do bi = itlo,ithi
                0434            do j = jmin,jmax
                0435             do i = imin,imax
                0436              junk = anom_sst(i,j,bi,bj)
f09238ab8f Gael*0437              junkweight = gencost_weight(i,j,bi,bj,kgen)*
f1e8c1d7ee Gael*0438      &          maskc(i,j,1,bi,bj)*msk_sst(i,j,bi,bj)
f4cfad5b95 Ou W*0439              objf_gencost(bi,bj,kgen) =
                0440      &          objf_gencost(bi,bj,kgen)+junk*junk*junkweight
f1e8c1d7ee Gael*0441              if (junkweight.GT.0.)
f4cfad5b95 Ou W*0442      &          num_gencost(bi,bj,kgen) =
                0443      &          num_gencost(bi,bj,kgen) + 1. _d 0
f1e8c1d7ee Gael*0444             enddo
                0445            enddo
                0446           enddo
                0447          enddo
                0448 
                0449       enddo
                0450 
f09238ab8f Gael*0451 c ------
                0452       endif !if (kgen.NE.0) then
                0453 
f1e8c1d7ee Gael*0454 #endif /* ifdef ALLOW_GENCOST_CONTRIBUTION */
                0455 #endif /* ifdef ALLOW_SMOOTH */
                0456 
f8e779c983 antn*0457       RETURN
                0458       END