PROGRAM extractrma11
    
c***************************************************************************************************
c This program reads in a series of RMA-11 results files and outputs the list of constituents
c at selected nodes.
c
c***************************************************************************************************

	IMPLICIT NONE

        CHARACTER*16 convtodate
        
        CHARACTER(len=60),dimension(:),allocatable :: rmafile
        INTEGER num_rma_files
        
        CHARACTER(len=100),dimension(:),allocatable :: output_name
        INTEGER,dimension(:),allocatable :: output_node
        INTEGER num_outputs
        
        INTEGER i,j,k,n,nqal
        
        INTEGER num_nodes, num_elems, num_qual
        
        CHARACTER(len=8),dimension(:),allocatable :: constit

        REAL,dimension(:,:),allocatable :: vel,tcon       
        REAL,dimension(:),allocatable :: wsel,wd

        INTEGER         IYRR
        REAL 		TETT
        
	CHARACTER(len=100) :: filename
	CHARACTER(len=40) :: out_prefix
	CHARACTER*1000	header


c -------------------------------------------------------------------------------------	
c       Open and Read the Config File
c--------------------------------------------------------------------------------------
	WRITE(*,'(A$)') 'Enter the config file name:'
	READ(*,'(A)') filename
        OPEN(UNIT=10,FILE=filename,STATUS='old')

        READ(10,'(I5)') num_rma_files
        allocate(rmafile(num_rma_files))
        do i=1,num_rma_files
                read(10,'(A)') rmafile(i)
        enddo
        
        READ(10,'(A)') out_prefix
        
        READ(10,'(I5)') num_outputs
        allocate(output_name(num_outputs))
        allocate(output_node(num_outputs))       
        
        do i=1,num_outputs
             read(10,'(I5,3X,A)') output_node(i), output_name(i)
             write(*,'(I5,A)') output_node(i),trim(output_name(i))
        enddo
 
        CLOSE(10)
        
c -------------------------------------------------------------------------------------
c   There is an assumption that all of the RMA files have come from the same mesh
c   and configuration.  You could test the header of each one, but let's assume the 
c   user knows what they are doing.
c
c   Open the first file and read the constituent list and number of nodes in the mesh
c   The allocate the arrays
c -------------------------------------------------------------------------------------

        OPEN(UNIT=20,FILE=rmafile(1),STATUS='old',FORM='binary')
        READ(20) header
        CLOSE(20)
        READ(header(41:70),'(3I10)') num_nodes,num_elems,num_qual
        allocate(constit(num_qual))
        DO i=1,num_qual
                read(header(302+8*(i-1):302+8*i),'(A)') constit(i)
        ENDDO
        allocate(vel(3,num_nodes))
        allocate(tcon(20,num_nodes))
        allocate(wd(num_nodes))
        allocate(wsel(num_nodes))
     
c -------------------------------------------------------------------------------------	
c       Open each CSV file and write the header
c--------------------------------------------------------------------------------------
        DO J=1,num_outputs
            write(filename,'(A,A,A)') trim(out_prefix),
     +                                trim(output_name(j)),'.csv'
            OPEN(UNIT=J+50,FILE=filename)
             
             write(J+50,'(A,$)'), 'Year,Time,Date,'
             DO I=1,num_qual
                 WRITE(J+50,('(A,$)')) constit(I),','
             ENDDO
             WRITE(J+50,('(A)')) ''   
        ENDDO

c -------------------------------------------------------------------------------------	
c       Loop through the RMA-11 files
c--------------------------------------------------------------------------------------
        
        DO i=1,num_rma_files
           write(*,'(A)') rmafile(i)
           OPEN(UNIT=12, FILE=rmafile(i),STATUS='OLD',FORM='BINARY')
           read(12) header
100        read(12,end=200) TETT,NQAL,N,IYRR,
     +          ((VEL(K,J),J=1,N),K=1,3), (wd(J),J=1,N),
     +          (wsel(J),J=1,N),
     +          ((TCON(K,J),J=1,N),K=1,NQAL-5)
     
           DO J=1,num_outputs              
                WRITE(J+50,('(I4,A,F7.2,A,A16,A$)')) 
     +                IYRR,',',TETT,',',convtodate(IYRR,TETT),','
                DO K=1,num_qual
                        WRITE(J+50,('(E14.7,A1,$)')) 
     +                              TCON(K,output_node(j)),','
                enddo
                WRITE(J+50,('(A)')) '' 
           ENDDO
	   GOTO 100
200	   CLOSE(12)
        ENDDO

c Close the files

        DO J=1,num_outputs
            CLOSE(UNIT=J+50)
        ENDDO

	END
	
c +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
c Subroutines
c +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

        CHARACTER*16 FUNCTION convtodate(IYRR,TETT)
        INTEGER IYRR
        REAL TETT
        INTEGER year,month,day,hour,minute,jday
        
        jday=jd(IYRR,1,1)+floor(TETT/24.0)
        CALL gdate(jday,year,month,day)
        hour=mod(floor(TETT),24)
        minute=int((TETT-floor(TETT))*60.0)
        
        write(convtodate,'(I2.2,A,I2.2,A,I4,A,I2.2,A,I2.2)')
     +                  day,'/',month,'/',year,' ',hour,':',minute
     
        RETURN
        END
        

        INTEGER FUNCTION JD (YEAR,MONTH,DAY)
C
C---    COMPUTES THE JULIAN DATE (JD) GIVEN A GREGORIAN CALENDAR
C       DATE (YEAR,MONTH,DAY).
C
        INTEGER YEAR,MONTH,DAY,I,J,K
C
        I= YEAR
        J= MONTH
        K= DAY
C
        JD= K-32075+1461*(I+4800+(J-14)/12)/4+367*(J-2-(J-14)/12*12)
     +            /12-3*((I+4900+(J-14)/12)/100)/4
C
        RETURN
        END
    
        SUBROUTINE GDATE (JD, YEAR,MONTH,DAY)
C
C---COMPUTES THE GREGORIAN CALENDAR DATE (YEAR,MONTH,DAY)
C   GIVEN THE JULIAN DATE (JD).
C
        INTEGER JD,YEAR,MONTH,DAY,I,J,K
C
        L= JD+68569
        N= 4*L/146097
        L= L-(146097*N+3)/4
        I= 4000*(L+1)/1461001
        L= L-1461*I/4+31
        J= 80*L/2447
        K= L-2447*J/80
        L= J/11
        J= J+2-12*L
        I= 100*(N-49)+I+L
C
        YEAR= I
        MONTH= J
        DAY= K
C
        RETURN
        END