C++*********************************************************************
C
C SSNRB.F
C              OPFILEC                             FEB  03 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 SSNRB(INC,Y1,WI,ILIST,NANG,NNNN,NSAM,NROW,LUN9)
C 
C IMAGE_PROCESSING_ROUTINE
C
C23456789012345678901234567890123456789012345678901234567890123456789012
C **********************************************************************

        SUBROUTINE SSNRB(INC,Y1,WI,ILIST,NANG,NNNN,NSAM,NROW,LUN9)
     
        INCLUDE 'CMBLOCK.INC'
        PARAMETER (NDLI=5)

        COMMON  /F_SPEC/  FINPAT,FINPIC,NLET
        DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: SUMCMP
        DOUBLE COMPLEX, ALLOCATABLE, DIMENSION(:,:) :: FA
        REAL, ALLOCATABLE, DIMENSION(:,:) :: B,RMSK

C       LR, SIGNAL AND FR ARE AUTOMATIC ARRAYS

        DIMENSION DLIST(NDLI), LR(INC),ILIST(NANG)
        DOUBLE PRECISION SIGNAL(INC),FR(INC)
        CHARACTER*80  FINPAT,FINPIC
 
        DATA  LUNI/99/

        ALLOCATE (RMSK(NSAM,NROW), STAT=IRTFLG)
        IF (IRTFLG.NE.0) THEN 
           CALL ERRT(46,'RF SN, RMSK',IER)
           RETURN
        ENDIF
  
        DO    J=1,NROW
          CALL REDLIN(LUN9,RMSK(1,J),NSAM,J)
        ENDDO
        CLOSE(LUN9)

        ALLOCATE (SUMCMP(NNNN/2,NROW), STAT=IRTFLG)
        IF (IRTFLG.NE.0) THEN 
           CALL ERRT(46,'RF SN, SUMCMP',IER)
           RETURN
        ENDIF

        ALLOCATE (FA(NNNN/2,NROW), STAT=IRTFLG)
        IF (IRTFLG.NE.0) THEN 
           CALL ERRT(46,'RF SN, FA',IER)
           RETURN
        ENDIF

        ALLOCATE (B(NNNN,NROW), STAT=IRTFLG)
        IF (IRTFLG.NE.0) THEN 
           CALL ERRT(46,'RF SN, B',IER)
           RETURN
        ENDIF


C       INITIALIZE THE SUMS
 
        SUMCMP=0.0D0
        FA=CMPLX(0.0,0.0)

        DO  K=1,NANG
C          READ ONE PROJECTION
           CALL  FILGET(FINPAT,FINPIC,NLET,ILIST(K),INTFLAG)
           MAXIM = 0
           CALL OPFILEC(0,.FALSE.,FINPIC,LUNI,'O',IFORM,LSAM,LROW,NSL,
     &             MAXIM,'DUMMY',.FALSE.,IRTFLG)

           IF (IRTFLG .NE. 0) THEN
               WRITE(NOUT,2032) FINPIC
2032           FORMAT(' FILE SKIPPED: ',A)
           ELSE
 
C             READ  IMAGE
              DO I=1,NROW
                 CALL  REDLIN(LUNI,B(1,I),NSAM,I)
              ENDDO
              CLOSE(LUNI)
c$omp         parallel do private(i,j)
              DO J=1,NROW
                 DO I=1,NSAM
                    B(I,J)=B(I,J)*RMSK(I,J)
                 ENDDO
              ENDDO

              INV=+1
              CALL  FMRS_2(B,NSAM,NROW,INV)
              IF(INV.NE.1)  THEN
                 CLOSE(LUNI)
                 CALL  ERRT(38,'FB',NE)
                 RETURN
              ENDIF

c$omp         parallel do private(i,j,ix)
              DO J=1,NROW
                 DO I=1,NNNN,2
                    IX=(I+1)/2
                    FA(IX,J)=FA(IX,J)+CMPLX(B(I,J),B(I+1,J))
                    SUMCMP(IX,J)=SUMCMP(IX,J)+
     &                 B(I,J)*DBLE(B(I,J))+B(I+1,J)*DBLE(B(I+1,J))
                 ENDDO
              ENDDO
           ENDIF
        ENDDO
      
        SIGNAL = 0.0D0
        FR     = 0.0D0
        LR     = 0

        DO 149 J=1,NROW
           JJ=J-1
           IF (JJ.GT.NROW/2)  JJ=JJ-NROW
           DO 149 I=1,NNNN/2
              IF (I.EQ.1.AND.JJ.LT.0) GOTO 149
              PII=0.5*SQRT((FLOAT(JJ)/FLOAT(NROW/2))**2+
     &        (FLOAT(I-1)/FLOAT(NSAM/2))**2)
              IF (PII.LE.0.5)  THEN
                 L=MIN0(MAX0(NINT(PII*Y1/WI)+1,1),INC)
                 LR(L)=LR(L)+1
                 SIGNAL(L)=SIGNAL(L)+ABS(FA(I,J))**2
                 FR(L)=FR(L)+SUMCMP(I,J)-ABS(FA(I,J))**2/NANG
              ENDIF
149     CONTINUE

C       SAVE RESULTS
        DO   L=1,INC
           IF(LR(L).NE.0)  THEN
              DLIST(1)=L
              DLIST(2)=FLOAT(L-1)/FLOAT(INC-1)*0.5
              DLIST(3)=
     &        AMAX1(0.0,SNGL((FLOAT(NANG-1)/NANG)*SIGNAL(L)/FR(L))-1.0)
              DLIST(4)=LR(L)
              N1=LR(L)
              N2=(NANG-1)*LR(L)
             DLIST(5)=SQRT(2.*((N2+N1-2)*(1+2*DLIST(3))+N1*DLIST(3)**2)
     &        /FLOAT(N1*(N2-4)))
              CALL  SAVD(LUN9,DLIST,NDLI,IRTFLG)
           ENDIF
        ENDDO

        CLOSE(LUN9)
        CALL SAVDC

        DEALLOCATE (SUMCMP) 
        DEALLOCATE (FA)
        DEALLOCATE (RMSK) 
        DEALLOCATE (B)
        END