C++********************************************************************* C C WPRO_N.F SPEEDED UP FEB 2000 ARDEAN LEITH C C ********************************************************************** C=* FROM: SPIDER - MODULAR IMAGE PROCESSING SYSTEM. AUTHOR: J.FRANK * C=* Copyright (C) 1985-2005 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 WPRO_N(B,NSAM,NROW,CUBE,LTB,IPCUBE,NN,PHI,THETA,PSI,RI) C C PURPOSE: COMPUTES PROJECTION(S) OF VOLUME ACCORDING TO EULERIAN ANGLES. C C IPCUBE: 1 - BEGINNING C 2 - END C 3 - IX C 4 - IY C 5 - IZ C C23456789012345678901234567890123456789012345678901234567890123456789012 C--********************************************************************* SUBROUTINE WPRO_N(B,NSAM,NROW,NSLICE,CUBE,LTB,IPCUBE, & NN,PHI,THETA,PSI,RI) REAL,DIMENSION(NSAM,NROW) :: B REAL,DIMENSION(NSAM,NROW,NSLICE) :: CUBE INTEGER,DIMENSION(5,NN) :: IPCUBE REAL,DIMENSION(9) :: DM DOUBLE PRECISION :: CPHI,SPHI,CTHE,STHE DOUBLE PRECISION :: CPSI,SPSI DOUBLE PRECISION :: QUADPI,DGR_TO_RAD DOUBLE PRECISION :: RAD_TO_DGR COMMON /PAR/ LDPX,LDPY,LDPZ PARAMETER (QUADPI = 3.1415926535897932384626) PARAMETER (DGR_TO_RAD = (QUADPI/180)) CPHI = DCOS(DBLE(PHI)*DGR_TO_RAD) SPHI = DSIN(DBLE(PHI)*DGR_TO_RAD) CTHE = DCOS(DBLE(THETA)*DGR_TO_RAD) STHE = DSIN(DBLE(THETA)*DGR_TO_RAD) CPSI = DCOS(DBLE(PSI)*DGR_TO_RAD) SPSI = DSIN(DBLE(PSI)*DGR_TO_RAD) DM(1) = CPHI*CTHE*CPSI-SPHI*SPSI DM(2) = SPHI*CTHE*CPSI+CPHI*SPSI DM(3) = -STHE*CPSI DM(4) = -CPHI*CTHE*SPSI-SPHI*CPSI DM(5) = -SPHI*CTHE*SPSI+CPHI*CPSI DM(6) = STHE*SPSI DM(7) = STHE*CPHI DM(8) = STHE*SPHI DM(9) = CTHE DM1 = DM(1) DM2 = DM(2) DM3 = DM(3) C ZERO THE WHOLE B ARRAY B = 0.0 DIM = MIN(NSAM,NROW,NSLICE) IF ((2*(RI+1) + 1) .LE. DIM) THEN C NO NEED TO CHECK IQX & IQY BOUNDARIES DO I=1,NN K=IPCUBE(4,I)+LDPY XB = IPCUBE(3,I)*DM(1) + IPCUBE(4,I)*DM(4) + & IPCUBE(5,I)*DM(7) + LDPX YB = IPCUBE(3,I)*DM(2) + IPCUBE(4,I)*DM(5) + & IPCUBE(5,I)*DM(8) + LDPY ZB = IPCUBE(3,I)*DM(3) + IPCUBE(4,I)*DM(6) + & IPCUBE(5,I)*DM(9) + LDPZ DO J=IPCUBE(1,I),IPCUBE(2,I) IOX = IFIX(XB) IOY = IFIX(YB) IOZ = IFIX(ZB) DX = XB-IOX DY = YB-IOY DZ = ZB-IOZ A1 = CUBE(IOX,IOY,IOZ) A2 = CUBE(IOX+1,IOY,IOZ) - A1 A3 = CUBE(IOX,IOY+1,IOZ) - A1 A4 = CUBE(IOX,IOY,IOZ+1) - A1 A5 = -A2 - CUBE(IOX,IOY+1,IOZ) + CUBE(IOX+1,IOY+1,IOZ) A61= - CUBE(IOX,IOY,IOZ+1) + CUBE(IOX+1,IOY,IOZ+1) A6 = -A2 + A61 A7 = -A3 - CUBE(IOX,IOY,IOZ+1) + CUBE(IOX,IOY+1,IOZ+1) A8 = -A5 - A61 - CUBE(IOX,IOY+1,IOZ+1) + & CUBE(IOX+1,IOY+1,IOZ+1) B(J,K) = B(J,K) + & A1 + DZ*(A4 + A6*DX + (A7 + A8*DX)*DY) + A3*DY + & DX*(A2 + A5*DY) XB = XB + DM1 YB = YB + DM2 ZB = ZB + DM3 ENDDO ENDDO ELSE C MUST CHECK IQX & IQY BOUNDARIES DO I=1,NN K=IPCUBE(4,I)+LDPY XB = IPCUBE(3,I)*DM(1) + IPCUBE(4,I)*DM(4) + & IPCUBE(5,I)*DM(7) + LDPX YB = IPCUBE(3,I)*DM(2) + IPCUBE(4,I)*DM(5) + & IPCUBE(5,I)*DM(8) + LDPY ZB = IPCUBE(3,I)*DM(3) + IPCUBE(4,I)*DM(6) + & IPCUBE(5,I)*DM(9) + LDPZ DO J=IPCUBE(1,I),IPCUBE(2,I) C CHECK FOR PIXELS OUT OF BOUNDS IOX = IFIX(XB) IF (IOX.LT.1 .OR. IOX.GE.NSAM) GOTO 2 IOY = IFIX(YB) IF (IOY.LT.1 .OR. IOY.GE.NROW) GOTO 2 IOZ = IFIX(ZB) IF (IOZ.LT.1 .OR. IOZ.GE.NSLICE) GOTO 2 IOX = IFIX(XB) IOY = IFIX(YB) IOZ = IFIX(ZB) DX = XB-IOX DY = YB-IOY DZ = ZB-IOZ A1 = CUBE(IOX,IOY,IOZ) A2 = CUBE(IOX+1,IOY,IOZ) - A1 A3 = CUBE(IOX,IOY+1,IOZ) - A1 A4 = CUBE(IOX,IOY,IOZ+1) - A1 A5 = -A2 - CUBE(IOX,IOY+1,IOZ) +CUBE(IOX+1,IOY+1,IOZ) A61= -CUBE(IOX,IOY,IOZ+1) + CUBE(IOX+1,IOY,IOZ+1) A6 = -A2 + A61 A7 = -A3 - CUBE(IOX,IOY,IOZ+1) +CUBE(IOX,IOY+1,IOZ+1) A8 = -A5 - A61 - CUBE(IOX,IOY+1,IOZ+1) + & CUBE(IOX+1,IOY+1,IOZ+1) B(J,K) = B(J,K)+ A1 + DZ*(A4 + A6*DX + & (A7 + A8*DX)*DY) + A3*DY + DX*(A2 + A5*DY) 2 XB = XB + DM1 YB = YB + DM2 ZB = ZB + DM3 ENDDO ENDDO ENDIF END #ifdef NEVER A1 = CUBE(IOX,IOY,IOZ) A2 = CUBE(IOX+1,IOY,IOZ) - CUBE(IOX,IOY,IOZ) A3 = CUBE(IOX,IOY+1,IOZ) - CUBE(IOX,IOY,IOZ) A4 = CUBE(IOX,IOY,IOZ+1) - CUBE(IOX,IOY,IOZ) A5 = CUBE(IOX,IOY,IOZ) - CUBE(IOX+1,IOY,IOZ) & - CUBE(IOX,IOY+1,IOZ)+ CUBE(IOX+1,IOY+1,IOZ) A6 = CUBE(IOX,IOY,IOZ) - CUBE(IOX+1,IOY,IOZ) & - CUBE(IOX,IOY,IOZ+1)+ CUBE(IOX+1,IOY,IOZ+1) A7 = CUBE(IOX,IOY,IOZ) - CUBE(IOX,IOY+1,IOZ) & - CUBE(IOX,IOY,IOZ+1)+ CUBE(IOX,IOY+1,IOZ+1) A8 = CUBE(IOX+1,IOY,IOZ) + CUBE(IOX,IOY+1,IOZ) & + CUBE(IOX,IOY,IOZ+1) & - CUBE(IOX,IOY,IOZ)- CUBE(IOX+1,IOY+1,IOZ) & - CUBE(IOX+1,IOY,IOZ+1) & - CUBE(IOX,IOY+1,IOZ+1) + CUBE(IOX+1,IOY+1,IOZ+1) B(J,K)=B(J,K)+ & A1 + DZ*(A4 + A6*DX + (A7 + A8*DX)*DY) + A3*DY & + DX*(A2 + A5*DY) #endif