C++********************************************************************* C SPEAK.F PARAMETERS FEB 2001 ARDEAN LEITH C FORMATTING MAY 2001 ARDEAN LEITH C NDAT/NOUT BUG MAY 2007 ARDEAN LEITH C ********************************************************************** C=* FROM: SPIDER - MODULAR IMAGE PROCESSING SYSTEM. AUTHOR: J.FRANK * C=* Copyright (C) 1985-2007 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 PURPOSE: SEARCHES FOR THE ML HIGHEST PEAKS IN THE (REAL) IMAGE C FILNAM AND PRINTS OUT POSITIONS AND VALUES OF THESE PEAKS. C C SPEAK(FILNAM,LUN,NSAM,NROW,MAXDIM,OPT,NDOC) c C PARAMETERS: C FILNAM FILE NAME C LUN LOGICAL UNIT NUMBER OF IMAGE C NSAM,NROW DIMENSIONS OF IMAGE C MAXDIM MAXIMUM BUFFER SPACE AVAILABLE C OPT OUTPUT OPTION C OPT=' ' DEFAULT:NO DOCUMENT OUTPUT C OPT='X' (I.E.,REGISTER LIST FOLLOWING):NO DOCUMENT OUTPUT C BUT OUTPUT OF POSITION & VALUE OF PEAK IN REGISTERS) C OPT='D' DOCUMENT OUTPUT:NUMBER,POSITION, AND VALUE C OF PEAK ARE WRITTEN INTO A DOCUMENT FILE C NDOC LOGICAL UNIT NUMBER FOR DOCUMENT FILE C C REGISTER POSITIONS 1= INTEGER X-SHIFT C 2= INTEGER Y-SHIFT C 3= ABSOLUTE PEAK HEIGHT C C 5= FLOATING PT. X-SHIFT C 6= FLOATING PT. Y-SHIFT C 7= VALUE OF EXTREMUM OF PARABOLOID C C23456789012345678901234567890123456789012345678901234567890123456789012 C--********************************************************************* SUBROUTINE SPEAK(FILNAM,LUN,NSAM,NROW,MAXDIM,OPT,NDOC,ML,NOR) INCLUDE 'CMBLOCK.INC' COMMON BUF(1) CHARACTER *(*) FILNAM CHARACTER OPT DIMENSION DLIST(5),RSQ(9) REAL, DIMENSION(ML) :: ALIST,RLIST INTEGER, DIMENSION(ML) :: ILIST,KLIST,KXLIST,IXLIST #ifdef USE_MPI include 'mpif.h' ICOMM = MPI_COMM_WORLD MPIERR = 0 CALL MPI_COMM_RANK(ICOMM, MYPID, MPIERR) #else MYPID = -1 #endif IF (MAXDIM .LT. 3*NSAM) THEN CALL ERRT(6,'SPEAK ',NE) RETURN ENDIF DO K = 1,ML ILIST(K) = 0 KLIST(K) = 0 RLIST(K) = 0.0 ALIST(K) = -HUGE(THRESH) ENDDO NTAB = 1 NXCTR = NSAM/2+1 NYCTR = NROW/2+1 IF (NOR .NE. 0) THEN CALL RDPRIS(NXCTR,NYCTR,NOT_USED, & 'ENTER ORIGIN COORDINATES',IRTFLG) IF (IRTFLG .NE. 0) RETURN CALL RDPRI1S(NTAB,NOT_USED, & 'ENTER PEAK NUMBER FOR RATIO',IRTFLG) IF (IRTFLG .NE. 0) RETURN IF (NTAB .LE. 0 .OR. NTAB .GT. ML) THEN CALL ERRT(25,'SPEAK',NE) RETURN ENDIF ENDIF 30 THRESH = -HUGE(THRESH) NMAX = 0 NSAM1 = NSAM-1 NROW1 = NROW -1 I1 = 0 I2 = 1 I3 = 2 CALL REDLIN(LUN,BUF,NSAM,1) CALL REDLIN(LUN,BUF(1+NSAM),NSAM,2) DO I = 3,NROW I1 = MOD(I1,3)+1 I2 = MOD(I2,3)+1 I3 = MOD(I3,3)+1 I1A = (I1-1)*NSAM I2A = (I2-1)*NSAM I3A = (I3-1)*NSAM CALL REDLIN(LUN,BUF(1+I3A),NSAM,I) DO 150 K = 2,NSAM1 A = BUF(K+I2A) C IGNORE PIXEL IF LOWER THAN LOWEST PIXEL ON PEAK LIST IF (A.LE.THRESH) GOTO 150 C IGNORE PIXEL IF LESS OR EQUAL TO ANY OF 8 NEIGHBORS IF (A.LE.BUF(K-1+I2A)) GOTO 150 IF (A.LE.BUF(K-1+I1A)) GOTO 150 IF (A.LE.BUF(K+I1A)) GOTO 150 IF (A.LE.BUF(K+1+I1A)) GOTO 150 IF (A.LE.BUF(K+1+I2A)) GOTO 150 IF (A.LE.BUF(K+1+I3A)) GOTO 150 IF (A.LE.BUF(K+I3A)) GOTO 150 IF (A.LE.BUF(K-1+I3A)) GOTO 150 NMAX = NMAX + 1 DO 100 L = 1,ML IF (A .LE. ALIST(L)) GOTO 100 IF (L .EQ. ML)GOTO 90 L1 = L+1 DO J = L1,ML JO = ML-J+L1 JN = ML-J+L1-1 ALIST(JO) = ALIST(JN) ILIST(JO) = ILIST(JN) KLIST(JO) = KLIST(JN) ENDDO 90 ALIST(L) = A ILIST(L) = I-1 KLIST(L) = K IF (NMAX .GT. ML) THRESH = ALIST(ML) GOTO 150 100 CONTINUE 150 CONTINUE ENDDO IF (NMAX .EQ. 0) THEN IF (MYPID .LE. 0) WRITE(NDAT,*) ' NO PEAK FOUND' IF (NDAT .NE. NOUT .AND. MYPID .LE. 0) & WRITE(NOUT,*) ' NO PEAK FOUND' CALL REG_SET_NSEL(1, 5,0.0, 0.0, 0.0, 0.0, 0.0,IRTFLG) CALL REG_SET_NSEL(6, 2,0.0, 0.0, 0.0, 0.0, 0.0,IRTFLG) RETURN ENDIF IF (MYPID .LE. 0) WRITE(NDAT,299) IF (NDAT .NE. NOUT .AND. MYPID .LE. 0) WRITE(NOUT,299) 299 FORMAT(' NO NSAM NROW X Y VALUE RATIO') MLIST = MIN0(ML,NMAX) DO L = 1,MLIST KXLIST(L) = KLIST(L) - NXCTR IXLIST(L) = ILIST(L) - NYCTR RLIST(L) = ALIST(L) / ALIST(NTAB) C IF (L .EQ. 1 .AND. NDAT .NE. NOUT) IF (OPT .EQ. 'D') THEN DLIST(1) = L DLIST(2) = KXLIST(L) DLIST(3) = IXLIST(L) DLIST(4) = ALIST(L) DLIST(5) = RLIST(L) CALL SAVD(NDOC,DLIST,5,IRTFLG) ENDIF IF (NDAT .NE. NOUT .AND. MYPID .LE. 0) & WRITE(NOUT,301)L,KLIST(L),ILIST(L),KXLIST(L),IXLIST(L), & ALIST(L),RLIST(L) 301 FORMAT(5I6,G16.7,2X,F8.5) IF (MYPID .LE. 0) & WRITE(NDAT,301) L,KLIST(L),ILIST(L),KXLIST(L), & IXLIST(L), ALIST(L),RLIST(L) ENDDO CALL SAVDC CLOSE(NDOC) IF (MYPID .LE. 0) WRITE(NDAT,302) IF (NDAT .NE. NOUT .AND. MYPID .LE. 0) WRITE(NOUT,302) 302 FORMAT(/,' LARGEST PEAK:') CALL REG_SET_NSEL(1,4,FLOAT(KXLIST(1)),FLOAT(IXLIST(1)), & ALIST(1),RLIST(1),0.0,IRTFLG) C 9/25/81 PARABOLIC FIT TO THE 3X3 NEIGHBORHOOD OF THE PEAK POINT C PROGRAM SECTION SENT BY M.VAN HEEL, MODIFIED FOR SPIDER USE. JF KL = KLIST(1) DO I=1,3 IROW = ILIST(1)+I-2 IF (IROW.LT.1) IROW=IROW+NROW IF (IROW.GT.NROW) IROW = IROW-NROW CALL REDLIN(LUN,BUF,NSAM,IROW) I1 = (I-1)*3 DO K=1,3 ISAM = KL+K-2 IF (ISAM.LT.1) ISAM=ISAM+NSAM IF (ISAM.GT.NSAM) ISAM=ISAM-NSAM RSQ(I1+K) = BUF(ISAM) ENDDO ENDDO CALL PARABL(RSQ,XSH,YSH,PEAKV) IF (NOUT .NE. NDAT .AND. MYPID .LE. 0) & WRITE(NOUT,351)XSH,YSH,PEAKV IF (MYPID .LE. 0) WRITE(NDAT,351)XSH,YSH,PEAKV 351 FORMAT(' + (',F5.2,', ',F5.2,') PEAK VALUE = ',G15.7) XT = XSH + KXLIST(1) YT = YSH + IXLIST(1) CALL REG_SET_NSEL(5,3,XT,YT,PEAKV, 0.0,0.0,IRTFLG) RETURN END