C++********************************************************************* C C ********************************************************************** C=* * C=* This file is part of: SPIDER - Modular Image Processing System. * C=* SPIDER System Authors: Joachim Frank & ArDean Leith * C=* Copyright 1985-2010 Health Research Inc., * C=* Riverview Center, 150 Broadway, Suite 560, Menands, NY 12204. * C=* Email: spider@wadsworth.org * C=* * C=* SPIDER 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=* SPIDER 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=* You should have received a copy of the GNU General Public License * C=* along with this program. If not, see * C=* * C ********************************************************************** C C COMP3D(LUN,LUN2) C C PURPOSE: C C PARAMETERS: C C23456789 123456789 123456789 123456789 123456789 123456789 123456789 12 C--********************************************************************* SUBROUTINE COMP3D(LUN,LUN2) IMPLICIT DOUBLE PRECISION (A-H,O-X) INCLUDE 'CMBLOCK.INC' INCLUDE 'CMLIMIT.INC' CHARACTER (LEN=MAXNAM) :: FIL2,FIL1 COMMON /COMMUN/ FIL2,FIL1 REAL, ALLOCATABLE, DIMENSION(:) :: AIMG INTEGER NVOX WRITE(NOUT,*) ' COMPARISON OF 2 3D ARRAYS:' MAXIM = 0 CALL OPFILEC(0,.TRUE.,FIL1,LUN,'O',IFORM1,NSAM,NROW,NSLICE, & MAXIM,'FIRST',.FALSE.,IRTFLG) IF (IRTFLG .NE. 0) GOTO 9900 MAXIM = 0 CALL OPFILEC(0,.TRUE.,FIL2,LUN2,'O',IFORM2,NSAM2,NROW2,NSLIC2, & MAXIM,'SECOND',.FALSE.,IRTFLG) IF (IRTFLG .NE. 0) GOTO 9900 IF(NSAM2.NE.NSAM.OR.NROW2.NE.NROW.OR.NSLIC2.NE.NSLICE) GOTO 9900 S1=0. S2=0. S3=0. S4=0. SUM1=0. SUM2=0. NVOX=FLOAT(NSAM)*FLOAT(NROW)*FLOAT(NSLICE) ALLOCATE (AIMG(2*NSAM), STAT=IRTFLG) IF (IRTFLG.NE.0) THEN CALL ERRT(46,'DR ERR, AIMG',IER) GOTO 9999 ENDIF DO ISL=1,NSLICE IOFF=(ISL-1)*NROW DO I = 1,NROW CALL REDLIN(LUN,AIMG,NSAM,I+IOFF) CALL REDLIN(LUN2,AIMG(NSAM+1),NSAM,I+IOFF) DO K = 1,NSAM B = AIMG(K) C=AIMG(NSAM+K) P1=B*C P2=C*C S1=S1+P1 S2=S2+P2 S3=S3+C S4=S4+B C IF((B.NE.0.OR.C.NE.0).AND.NSAM.LE.16) C $ WRITE(NOUT,103)ISL,I,K,B,C,P1,P2,S1,S2,S3,S4 C103 FORMAT(' SLICE: ',I3,' , ROW: ',I3,' , COL: ',I3,/ C $ ' VOXEL 1: ',F16.4,' , VOXEL 2: ',F16.4, / C $ ' AD1: ',F16.4,' AD2: ',F16.4, / C $ ' S1: ',F16.4,' S2: ',F16.4, / C $ ' S3: ',F16.4,' S4: ',F16.4) ENDDO ENDDO ENDDO A=(S1-S3*S4/NVOX)/(S2-S3*S3/NVOX) CONST=(A*S3-S4)/NVOX BMEAN=S4/NVOX SUM1=0 SUM2=0 DNOM=0 DNUM=0 DO ISL=1,NSLICE IOFF=(ISL-1)*NROW DO I = 1,NROW CALL REDLIN(LUN,AIMG,NSAM,I+IOFF) CALL REDLIN(LUN2,AIMG(NSAM+1),NSAM,I+IOFF) DO K = 1,NSAM B = AIMG(K) C=AIMG(NSAM+K) BD1=DABS(A*C-B-CONST) BD2=DABS(B-BMEAN) SUM1=SUM1+BD1 SUM2=SUM2+BD2 BD1SQ=BD1*BD1 DNOM=DNOM+(B-BMEAN)*(B-BMEAN) DNUM=DNUM+BD1SQ C IF((B.NE.0.OR.C.NE.0).AND.NSAM.LE.16) C $ WRITE(NOUT,104)ISL,I,K,B,C,BD1,BD2,SUM1,SUM2 C104 FORMAT(' SLICE: ',I3,' , ROW: ',I3,' , COL: ',I3,/ C $ ' VOXEL 1: ',F16.4,' , VOXEL 2: ',F16.4, / C $ ' BD1: ',F16.4,' BD2: ',F16.4, / C $ ' SUM1: ',F16.4,' SUM2: ',F16.4) ENDDO ENDDO ENDDO DIF = (SUM1)/SUM2 DISCR = DSQRT(DNUM/DNOM) WRITE(NOUT,100) A,CONST,DIF,DISCR 100 FORMAT(' SCALING OF SECOND 3D: VOXEL_NEW=VOXEL_OLD*',G16.6, & ' - ',G16.6/ & ' MEAN RELATIVE ERROR: ',G10.5/ & ' DISCREPANCY: ',G15.5) DEALLOCATE (AIMG) GOTO 9999 9900 CALL ERRT(IER,'COMP3D',NE) 9999 CLOSE (LUN) CLOSE (LUN2) END