C++********************************************************************* C C MRQLI_SS.F USED CMLIMIT AUG 00 ARDEAN LEITH C ADDED REF_CIRC FILE APR 01 ARDEAN LEITH C NORMASS -> NORMAS OCT 01 ARDEAN LEITH C PROMPTS JAN 02 ARDEAN LEITH C OPFILEC FEB 03 ARDEAN LEITH C APMASTER REWRITE AUG 03 ARDEAN LEITH C 'OR MQ' SUPPORT SEP 03 ARDEAN LEITH C AP_OUT USAGE FEB 04 ARDEAN LEITH C 'MQ' NOT READ EXP. ANGLES APR 04 ARDEAN LEITH C BAD PEAK IF INTERP ON EDGE AUG 04 ARDEAN LEITH C CROSRNG_E SPEEDS UP AUG 04 ARDEAN LEITH C AP_END CALL HAS PARLIST OCT 04 ARDEAN LEITH C LIMITRANGE BUG OCT 04 ARDEAN LEITH C PEAKV = 1 JAN 05 ARDEAN LEITH C DISCARD MIRROR ... JUN 06 ARDEAN LEITH C AP_STAT CALL JAN 07 ARDEAN LEITH C USE FFTW3 IN APRINGS MAR 08 ARDEAN LEITH C APRINGS_ONE_NEW APR 08 ARDEAN LEITH C C ********************************************************************** C=* FROM: SPIDER - MODULAR IMAGE PROCESSING SYSTEM. AUTHOR: J.FRANK * C=* Copyright (C) 1985-2008 Health Research Inc. * C=* * C=* HEALTH RESEARCH INCORPORATED (HRI), * C=* ONE UNIVERSITY PLACE, RENSSELAER, NY 12144-3455. * C=* * C=* Email: spider@wadsworth.org * C=* * C=* This program is free software; you can redistribute it and/or * C=* modify it under the terms of the GNU General Public License as * C=* published by the Free Software Foundation; either version 2 of the * C=* License, or (at your option) any later version. * C=* * C=* This program is distributed in the hope that it will be useful, * C=* but WITHOUT ANY WARRANTY; without even the implied warranty of * C=* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU * C=* General Public License for more details. * C=* * C=* You should have received a copy of the GNU General Public License * C=* along with this program; if not, write to the * C=* Free Software Foundation, Inc., * C=* 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. * C=* * C ********************************************************************** C C MRQLI_SS.F C C PURPOSE: FIND ROTATIONAL AND SHIFT PARAMETERS TO ALIGN A SERIES OF C REFERENCE IMAGES WITH SAMPLE IMAGES C C VERSION FOR MP AND A SMALL NUMBER OF IMAGES TO BE ALIGNED C C23456789012345678901234567890123456789012345678901234567890123456789012 C--********************************************************************* SUBROUTINE MRQLI_SS(IREFLIST,NUMREF,IEXPLIST,NUMEXP, & LSAM,LROW,NR,LENTT,ISHRANGE,ISTEP, & NRING,LCIRC,NUMR,CIRCREF,CIRCREF_IN_CORE, & MODE,REFANGDOC,EXPANGDOC,SCRFILE,FFTW_PLANS, & REFPAT,EXPPAT,RANGE,CKMIRROR,CTYPE,LUNDOC) INCLUDE 'CMBLOCK.INC' INCLUDE 'CMLIMIT.INC' INTEGER, DIMENSION(NUMREF) :: IREFLIST INTEGER, DIMENSION(NUMEXP) :: IEXPLIST INTEGER, DIMENSION(3,NRING) :: NUMR REAL, DIMENSION(LCIRC,NUMREF) :: CIRCREF LOGICAL :: CIRCREF_IN_CORE,CKMIRROR LOGICAL :: MIRRORNEW LOGICAL :: GOTREFANG LOGICAL :: LIMITRANGE CHARACTER(LEN=1) :: MODE,NULL CHARACTER (LEN=*) :: REFANGDOC,EXPANGDOC CHARACTER (LEN=*) :: SCRFILE CHARACTER (LEN=*) :: REFPAT,EXPPAT CHARACTER (LEN=*) :: CTYPE DOUBLE PRECISION :: FITP(-1:1,-1:1) DOUBLE PRECISION :: CCROTD,PEAK,CCROTD_INTERP INTEGER *8 :: FFTW_PLANS(*) #ifdef __ia64 C INTEL PARALLEL COMPILER BUGGY MUST ALLOCATE THESE ARRAYS IF || C DOES NOT HANDLE NEGATIVE SUBSCRIPTS CORRECTLY DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: TOTA,TMTA REAL, ALLOCATABLE, DIMENSION(:,:,:) :: TOT,TMT DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: FIT REAL, ALLOCATABLE, DIMENSION(:,:) :: ROTMP #else C AUTOMATIC ARRAYS DOUBLE PRECISION TOTA(NUMREF,-ISHRANGE:ISHRANGE, & -ISHRANGE:ISHRANGE) DOUBLE PRECISION TMTA(NUMREF,-ISHRANGE:ISHRANGE, & -ISHRANGE:ISHRANGE) DIMENSION TOT(NUMREF, -ISHRANGE:ISHRANGE, & -ISHRANGE:ISHRANGE), & TMT(NUMREF, -ISHRANGE:ISHRANGE, & -ISHRANGE:ISHRANGE) DOUBLE PRECISION FIT( -ISTEP:ISTEP,-ISTEP:ISTEP) DIMENSION ROTMP(-ISTEP:ISTEP,-ISTEP:ISTEP) #endif DOUBLE PRECISION, DIMENSION(LENTT) :: TT INTEGER, DIMENSION(NUMREF) :: LCG REAL, DIMENSION(3) :: ANGOUT REAL, DIMENSION(3) :: EXPDIR LOGICAL :: MIRRORED C ALLOCATED ARRAYS REAL, ALLOCATABLE, DIMENSION(:,:,:) :: CIRCEXP REAL, ALLOCATABLE, DIMENSION(:,:) :: A REAL, ALLOCATABLE,DIMENSION(:,:) :: REFDIR REAL, ALLOCATABLE,DIMENSION(:,:) :: ANGREF,ANGEXP INTEGER, PARAMETER :: NLISTMAX = 15 REAL, DIMENSION(NLISTMAX) :: PARLIST PARAMETER (QUADPI = 3.1415926535897932384626) PARAMETER (DGR_TO_RAD = (QUADPI/180)) DATA INPIC/77/,INANG/78/,LUNRING/50/ NULL = CHAR(0) C INITIALIZE CCROT STATISTICS COUNTERS ANGDIFTHR = 0.0 IBIGANGDIF = -1 ANGDIFAVG = 0.0 CCROTAVG = 0.0 IMPROVCCROT = 0 CCROTIMPROV = 0.0 IWORSECCROT = 0 CCROTWORSE = 0.0 C FLAG FOR RESTRICTED PROJECTION RANGE LIMITRANGE = RANGE .GT. 0 .AND. RANGE .LT. 360 MAXRIN = NUMR(3,NRING) #ifdef SP_LIBFFTW3 MAXRIN = NUMR(3,NRING) -2 ! ACTUAL LENGTH OF LONGEST RING #endif RANGE = COS(RANGE*DGR_TO_RAD) WR = 0.0 ! DUMMY VALUE FLAG FOR APRINGS CALL C FIND DIVAS, NUMTH, NSAM, & NROW CALL APMASTER_1(MODE,DIVAS,NR,NUMTH,LSAM,LROW,NSAM,NROW, & TT,LENTT) C READ REFERENCE IMAGES INTO REFERENCE RINGS (CIRCREF) ARRAY CALL APRINGS_NEW(IREFLIST,NUMREF, LSAM,LROW, & NRING,LCIRC,NUMR,MODE,FFTW_PLANS, & REFPAT,INPIC,CIRCREF,CIRCREF_IN_CORE, & LUNRING,SCRFILE,IRTFLG) NSAM = LSAM NROW = LROW NSIS = MAX(ISHRANGE/ISTEP, ISTEP) ALLOCATE(CIRCEXP(LCIRC,-NSIS:NSIS,-NSIS:NSIS), & A(NSAM,NROW), & STAT=IRTFLG) IF (IRTFLG .NE. 0) THEN NSIST = 2 * NSIS + 1 MWANT = LCIRC*NSIST*NSIST + NSAM*NROW CALL ERRT(46,' CIRCEXP,....',MWANT) GOTO 9999 ENDIF NUMREFLCG = NUMREF IEND = NUMREF NGOTPAR = 0 IF (LIMITRANGE .OR. CTYPE(1:2) .EQ. 'SH') THEN C RESTRICTED ANGULAR SEARCH FOR 'RQ' OR 'SH' ALLOCATE(ANGREF(3,NUMREF), STAT=IRTFLG) IF (IRTFLG .NE. 0) THEN MWANT = 3*NUMREF CALL ERRT(46,'ANGREF',MWANT) RETURN ENDIF C LOAD REF. PROJ. ANGLES () FROM DOC. FILE (REFANGDOC) OR C IMAGE FILE (REFPAT) HEAD INTO ANGREF CALL AP_GETANGA(IREFLIST,NUMREF,0,REFANGDOC,REFPAT, & INPIC,INANG,3,ANGREF,NGOTREF,IRTFLG) IF (IRTFLG .NE. 0) GOTO 9999 ENDIF IF (CTYPE(1:2) .NE. 'MQ' .AND. & EXPANGDOC .NE. CHAR(0)) THEN C READ EXP. ANGLES INTO ANGEXP IF (REFANGDOC .NE. NULL) THEN ALLOCATE(REFDIR(3,NUMREF),STAT=IRTFLG) IF (IRTFLG .NE. 0) THEN MWANT = 3*NUMREF CALL ERRT(46,'REFDIR',MWANT) RETURN ENDIF C CONVERT REF. ANGLES TO UNITARY DIRECTIONAL VECTORS (REFDIR). CALL AP_GETSATA(ANGREF,REFDIR,3,NUMREF,IRTFLG) ENDIF C LOAD EXP. PROJ. ANGLES & ALIGNMENT PARAMETERS (ANGEXP) C FROM DOC. FILE (EXPANGDOC) OR IMAGE FILE (REFPAT) HEAD ALLOCATE(ANGEXP(7,NUMEXP), STAT=IRTFLG) IF (IRTFLG .NE. 0) THEN MWANT = 7*NUMEXP CALL ERRT(46,'ANGEXP',MWANT) RETURN ENDIF C THIS RETURNS NGOTPAR CALL AP_GETANGA(IEXPLIST,NUMEXP,0,EXPANGDOC,EXPPAT, & INPIC,INANG,7,ANGEXP,NGOTPAR,IRTFLG) IF (IRTFLG .NE. 0) GOTO 9999 GOTREFANG = .TRUE. ENDIF C LOOP OVER EXP. IMAGES TO BE ALIGNED DO IEXP=1,NUMEXP IMGEXP = IEXPLIST(IEXP) IF (CTYPE(1:2) .EQ. 'RQ' .OR. & (LIMITRANGE .AND. EXPANGDOC .NE. NULL)) THEN C CONVERT EXP. ANGLE TO UNITARY DIRECTIONAL VECTORS (EXPDIR). CALL AP_GETSATA(ANGEXP(1,IEXP),EXPDIR,7,1,IRTFLG) C MAKE LIST OF NEARBY REFERENCE IMAGES NUMREFLCG = 0 DO IMIR=1,NUMREF C LOOP OVER REFERENCE IMAGES C DT NEAR 1.0 = NOT-MIRRORED, DT NEAR -1.0 = MIRRORED DT = (EXPDIR(1) * REFDIR(1,IMIR) + & EXPDIR(2) * REFDIR(2,IMIR) + & EXPDIR(3) * REFDIR(3,IMIR)) DTABS = ABS(DT) IF (DTABS .GE. RANGE) THEN C MIRRORED OR NON-MIRRORED IS WITHIN RANGE IF (CKMIRROR .OR. DT .GT. 0) THEN C DO NOT DISCARD IF NOT MIRRORED OR WANT MIRRORED NUMREFLCG = NUMREFLCG + 1 LCG(NUMREFLCG) = IMIR IF (DT .LT. 0) LCG(NUMREFLCG) = -IMIR ENDIF IEND = NUMREFLCG ENDIF ENDDO IF (NUMREFLCG .LE. 0) THEN C REPORT THAT THERE ARE NO NEARBY REFERENCE IMAGES IMGREF = 0 PEAKV = 0.0 CCROT = 0.0 XSHNEW = 0.0 YSHNEW = 0.0 MIRRORNEW = .FALSE. CALL AP_END(IEXP,IMGEXP,IMGREF, & ANGREF,REFDIR, & ANGEXP(1,IEXP),EXPDIR,ISHRANGE, & GOTREFANG, NGOTPAR,LSAM,LROW,CCROT,PEAKV, & RANGNEW,XSHNEW,YSHNEW,MIRRORNEW,EXPPAT,REFPAT, & NUMREFLCG, CTYPE, XBUF,LUNDOC,PARLIST) CYCLE ENDIF ENDIF C LOAD EXP. IMAGE DATA FOR THIS IMAGE IMITT = IEXP IF (CTYPE(1:2) .EQ. 'OR') IMITT = 0 CALL AP_GETDATS(IEXPLIST,NUMEXP,LSAM,LROW, & NUMTH,EXPPAT,INPIC,IMITT,IMITT, & A, IRTFLG) IF (IRTFLG .NE. 0) GOTO 9999 #ifdef __ia64 ALLOCATE(TOTA(NUMREF,-ISHRANGE:ISHRANGE, -ISHRANGE:ISHRANGE), & TMTA(NUMREF,-ISHRANGE:ISHRANGE, -ISHRANGE:ISHRANGE), & TOT(NUMREF, -ISHRANGE:ISHRANGE, -ISHRANGE:ISHRANGE), & TMT(NUMREF, -ISHRANGE:ISHRANGE, -ISHRANGE:ISHRANGE), & FIT( -ISTEP:ISTEP, -ISTEP:ISTEP), & ROTMP(-ISTEP:ISTEP, -ISTEP:ISTEP), & STAT=IRTFLG) IF (IRTFLG .NE. 0) THEN MWANT1 = 4 * (NUMREF * (ISHRANGE*2+1) * (ISHRANGE*2+1)) + & 2 * (ISTEP*2+1) CALL ERRT(46,' TOTA,....',MWANT1) GOTO 9999 ENDIF #endif C GO THROUGH REFERENCE IMAGES AND CENTERS c$omp parallel do private(jt,it,irr,ir,cnr2,cns2,mirrored) DO JT=-ISHRANGE,ISHRANGE,ISTEP CNR2 = NROW/2+1+JT DO IT=-ISHRANGE,ISHRANGE,ISTEP CNS2 = NSAM/2+1+IT C NORMALIZE IMAGE VALUES UNDER THE MASK OVER VARIANCE RANGE C INTERPOLATE TO POLAR COORDINATES, CREATE FOURIER OF: CIRCEXP C NO WEIGHTING OF RINGS CALL APRINGS_ONE_NEW(NSAM,NROW, CNS2,CNR2, A,.FALSE., & MODE,NUMR,NRING,LCIRC,WR,FFTW_PLANS, & CIRCEXP(1,IT/ISTEP,JT/ISTEP),IRTFLG) DO IRR=1,IEND C LOOP OVER REFERENCE IMAGES IR = IRR IF (LIMITRANGE) IR = ABS(LCG(IRR)) IF ((CKMIRROR .AND. LIMITRANGE) .OR. & (.NOT. CKMIRROR)) THEN IF (.NOT. CKMIRROR) MIRRORED = .FALSE. IF (LIMITRANGE) MIRRORED = (LCG(IRR) .LT. 0) C CHECK EITHER MIRRORED OR NON-MIRRORED POSITION CALL CROSRNG_EP_NEW(CIRCREF(1,IR), & CIRCEXP(1,IT/ISTEP,JT/ISTEP), & LCIRC,NRING,MAXRIN,NUMR, & TOTA(IR,IT,JT),TOT(IR,IT,JT), & TT,MIRRORED, FFTW_PLANS(1)) ELSE C CHECK BOTH NON-MIRRORED & MIRRORED POSITIONS CALL CROSRNG_MSP_NEW(CIRCREF(1,IR), & CIRCEXP(1,IT/ISTEP,JT/ISTEP), & LCIRC,NRING,MAXRIN,NUMR, & TOTA(IR,IT,JT),TOT(IR,IT,JT), & TMTA(IR,IT,JT),TMT(IR,IT,JT), & TT,FFTW_PLANS(1)) ENDIF ENDDO ! END OF: DO IRR=1,IEND ENDDO ! END OF: DO IT=-ISHRANGE,ISHRANGE,ISTEP ENDDO ! END OF: DO JT=-ISHRANGE,ISHRANGE,ISTEP C END OF THE OMP PARALLEL SECTION C LOCATE BEST CCROT MATCH FROM ALL THE VALUES ACCUMULATED ABOVE CCROTD = -1.0D23 DO JT=-ISHRANGE,ISHRANGE,ISTEP DO IT=-ISHRANGE,ISHRANGE,ISTEP DO IRR=1,IEND C LOOP OVER REFERENCE IMAGES IR = IRR IF (LIMITRANGE) IR = ABS(LCG(IRR)) IF (TOTA(IR,IT,JT) .GE. CCROTD) THEN C BETTER POSITION (MAY BE MIRRORED IF LIMITRANGE) CCROTD = TOTA(IR,IT,JT) IREF = IR ISX = IT ISY = JT RANGNEW = ANG_N(TOT(IR,IT,JT),MODE,MAXRIN) MIRRORNEW = (LIMITRANGE .AND. (LCG(IRR) .LT. 0)) ENDIF IF (CKMIRROR .AND. .NOT. LIMITRANGE) THEN C HAVE TO COMPARE WITH MIRRORED POSITION IF (TMTA(IR,IT,JT) .GE. CCROTD) THEN C BETTER MIRRORED POSITION IREF = IR CCROTD = TMTA(IR,IT,JT) ISX = IT ISY = JT RANGNEW = ANG_N(TMT(IR,IT,JT),MODE,MAXRIN) MIRRORNEW = .TRUE. ENDIF ENDIF ENDDO ENDDO ENDDO C write(6,*) 'iref: ',iref, ' ccrotd: ',ccrotd SX = ISX SY = ISY CCROT = CCROTD IMGREF = IREFLIST(IREF) #ifdef DEBUG write(6,921) IMGREF,isx,isy,ccrotd,rangnew 921 format(' 1 ',i5,' (',i3,',',i3,'): ',f14.4,' ',2f8.2,f6.1) #endif C CHECK LOCATIONS WITHIN ISHRANGE AROUND MAX ------------------ IF (ABS(ISX) .NE. ISHRANGE .AND. & ABS(ISY) .NE. ISHRANGE) THEN C NOT ON BOUNDARY, HAVE TO FIND NEIGHBOURING VALUES FIT(0,0) = CCROTD ROTMP(0,0) = RANGNEW c$omp parallel do private(jt,it,cnr2,cns2) DO JT=-ISTEP,ISTEP DO IT=-ISTEP,ISTEP CNR2 = NROW / 2 + 1 + JT + ISY IF (IT.NE.0 .OR. JT.NE.0) THEN CNS2 = NSAM / 2 + 1 + IT + ISX C NORMALIZE IMAGE VALUES UNDER THE MASK OVER VARIANCE RANGE C INTERPOLATE INTO POLAR COORDINATES C CREATES FOURIER OF: CIRCEXP CALL APRINGS_ONE_NEW(NSAM,NROW, CNS2,CNR2, & A,.FALSE., & MODE,NUMR,NRING,LCIRC,WR,FFTW_PLANS, & CIRCEXP(1,IT,JT),IRTFLG) CALL CROSRNG_EP_NEW(CIRCREF(1,IREF), & CIRCEXP(1,IT,JT), & LCIRC,NRING,MAXRIN,NUMR, & FIT(IT,JT),ROTMP(IT,JT), & TT, MIRRORNEW,FFTW_PLANS(1)) C RECORD BEST ANGLE IN ROTMP ROTMP(IT,JT) = ANG_N(ROTMP(IT,JT),MODE,MAXRIN) ENDIF ENDDO ! END OF: DO IT=-ISTEP,ISTEP ENDDO ! END OF: DO JT=-ISTEP,ISTEP C FIND THE MAXIMUM CC ANGLE WITHIN +/-ISTEP C MAXIMUM CANNOT BE ON THE EDGE, I.E., IT,JT/=ISTEP AFIT = FIT(0,0) JTMA = 0 ITMA = 0 RANGNEWT = ROTMP(0,0) IF (ISTEP .GT. 1) THEN DO JT=-ISTEP+1,ISTEP-1 DO IT=-ISTEP+1,ISTEP-1 IF (FIT(IT,JT) .GT. AFIT) THEN AFIT = FIT(IT,JT) RANGNEWT = ROTMP(IT,JT) !compiler bug on OPT64 ITMA = IT JTMA = JT ENDIF ENDDO ENDDO ENDIF ! END OF: IF (ISTEP .GT. 1) C write(6,*) ((fit(i,j),i=-1,1),j=-1,1) C TEMP VARIABLE OVERCOMES COMPILER BUG ON OPT 64 PGI 6.0 RANGNEW = RANGNEWT C COPY VALUES AROUND THE PEAK. DO JT=-1,1 DO IT=-1,1 FITP(IT,JT) = FIT(ITMA+IT,JTMA+JT) c write(6,910) it,jt,fitp(it,jt) 910 format(' fitp(',i5,',',i5,') : ',f14.4) ENDDO ENDDO C UPDATE INTEGRAL LOCATION OF THE PEAK CCROTD = AFIT ISX = ISX + ITMA ISY = ISY + JTMA SX = ISX SY = ISY #ifdef DEBUG cnr2 = nrow / 2 + 1 + isy cns2 = nsam / 2 + 1 + isx write(6,905)imgref,isx,isy,ccrotd,rangnew,cns2,cnr2,sx,sy 905 format(' 2 ',i5,' (',i3,',',i3,'): ',f12.4,' ',5f8.2) #endif C SUB-PIXEL INTERPOLATION ------------------------------ C FIND PEAK BY FITTING PARABOLA TO 3x3 NEIGHBORHOOD CALL PARABLD(FITP,SSX,SSY,PEAK) IF (ABS(SSX) .LT. 1.0 .AND. ABS(SSY) .LT. 1.0) THEN C INTERPOLATED LOCATION IS NOT ON BOUNDARY CNS2 = NSAM/2+1 + SX + SSX CNR2 = NROW/2+1 + SY + SSY C NORMALIZE IMAGE VALUES UNDER THE MASK OVER VARIANCE RANGE C INTERPOLATE INTO POLAR COORDINATES C CREATE FOURIER OF: CIRCEXP CALL APRINGS_ONE_NEW(NSAM,NROW, CNS2,CNR2, A,.FALSE., & MODE,NUMR,NRING,LCIRC,WR,FFTW_PLANS, & CIRCEXP,IRTFLG) CALL CROSRNG_EP_NEW(CIRCREF(1,IREF),CIRCEXP,LCIRC, & NRING,MAXRIN,NUMR,CCROTD_INTERP, & RANGNEW_INTERP,TT,MIRRORNEW,FFTW_PLANS(1)) #ifdef DEBUG rt1 = ang_n(rangnew_interp,mode,maxrin) write(6,904) imgref,isx,isy, ccrotd_interp,rt1, & cns2,cnr2, sx+ssx,sy+ssy 904 format(' 3 ',i5,' (',i3,',',i3,'): ',f12.4,' ',5f8.2) #endif IF (CCROTD_INTERP .GT. CCROTD) THEN C USE SUB-PIXEL LOCATION CCROTD = CCROTD_INTERP RANGNEW = ANG_N(RANGNEW_INTERP,MODE,MAXRIN) SX = SX + SSX SY = SY + SSY ENDIF ENDIF ! END OF: IF (ABS(SX) .LT. 1.0 ..... ENDIF ! END OF: IF ( ABS(ISX) .NE. ISHRANGE ...... CCROT = CCROTD SX = -SX SY = -SY C HAVE TO CHANGE ORDER OF SHIFT & ROTATION. C IN THIS PROGRAM IMAGE IS SHIFTED FIRST, ROTATED SECOND. C IN 'RT SQ' IT IS ROTATION FIRST, SHIFT SECOND. C THIS CODE CORRESPONDS TO 'SA P'. CO = COS(RANGNEW * DGR_TO_RAD) SO = -SIN(RANGNEW * DGR_TO_RAD) XSHNEW = SX*CO - SY*SO YSHNEW = SX*SO + SY*CO PEAKV = 1.0 C ALMOST ZERO IS LIKELY TO BE ZERO IF (ABS(XSHNEW) .LT. 0.04) XSHNEW = 0.0 IF (ABS(YSHNEW) .LT. 0.04) YSHNEW = 0.0 IF (ABS(RANGNEW) .LT. 0.04) RANGNEW = 0.0 CALL AP_END(IEXP,IMGEXP,IMGREF, & ANGREF,REFDIR(1,IREF), & ANGEXP(1,IEXP),EXPDIR,ISHRANGE, & GOTREFANG, NGOTPAR,LSAM,LROW,CCROT,PEAKV, & RANGNEW,XSHNEW,YSHNEW,MIRRORNEW,EXPPAT,REFPAT, & NUMREFLCG, CTYPE, XDUM,LUNDOC,PARLIST) CCROTAVG = CCROTAVG + CCROT IF (NGOTPAR .GE. 8) THEN C COMPILE CCROT CHANGE STATISTICS ANGDIF = PARLIST(10) IF (ANGDIF .GT. ANGDIFTHR)IBIGANGDIF = IBIGANGDIF + 1 CCROTLAS = ANGEXP(8,IEXP) ANGDIFAVG = ANGDIFAVG + PARLIST(10) IF (CCROT .GE. CCROTLAS) THEN IMPROVCCROT = IMPROVCCROT + 1 CCROTIMPROV = CCROTIMPROV + CCROT ELSE IWORSECCROT = IWORSECCROT + 1 CCROTWORSE = CCROTWORSE + CCROT ENDIF ENDIF ! END OF: IF (NGOTPAR .GE. 8) #ifdef DEBUG cns2 = nsam / 2 + 1 - sx cnr2 = nrow / 2 + 1 - sy write(6,906) imgref,isx,isy, ccrotd,rangnew, & cns2,cnr2, xshsum,yshsum 906 format(' 4 ',i5,' (',i3,',',i3,'): ',f12.4,' ',5f8.2) write(6,*) ' ------------------------------------- ' write(6,*) ' ' #endif ENDDO CALL REG_GET_USED(NSEL_USED) IF (NSEL_USED .GT. 0 .AND. CTYPE(1:2) .EQ. 'OR') THEN C OUTPUT TO REGISTER NOT TO DOC FILE (FOR 'OR') DMR = 0 IF (MIRRORNEW) DMR = 1 CALL REG_SET_NSEL(1,5,RANGNEW,XSHNEW,YSHNEW, & DMR,CCROT,IRTFLG) ENDIF IF (LUNDOC .GT. 0 .AND. NUMEXP .GT. 1) THEN C SAVE CCROT & ANGULAR DISPLACEMENT STATISTICS CALL AP_STAT(NUMEXP,ANGDIFTHR,IBIGANGDIF, & ANGDIFAVG, CCROTAVG, & IMPROVCCROT,CCROTIMPROV, & IWORSECCROT,CCROTWORSE,LUNDOC) ENDIF 9999 IF (ALLOCATED(CIRCEXP)) DEALLOCATE(CIRCEXP) IF (ALLOCATED(A)) DEALLOCATE(A) IF (ALLOCATED(ANGREF)) DEALLOCATE(ANGREF) IF (ALLOCATED(ANGEXP)) DEALLOCATE(ANGEXP) IF (ALLOCATED(REFDIR)) DEALLOCATE(REFDIR) END