CHARACTER*40 INFILE,FILENAME
	LOGICAL EFLG
	real*8 M1(3,4),M2(3,4),M3(3,4),M(3,4),ACC(3,4)
	REAL*8 P1(3),P2(3),T(3)
	character*40 string
	INTEGER*4 MM,NN,I,J

c	count=1
c 	read first operator into m1
	write(6,*) 'x2= T2 + [rot2] {T1+[rot1]X1}'
5	WRITE(6,*)'ENTER NAME OF right hand (rot1,t1) rot-trans op (.o):'
	READ(5,1200) INFILE
1200	FORMAT(A)
	IF (INFILE(:5).EQ.'     ') STOP
	  INQUIRE(FILE=INFILE, NAME=FILENAME,EXIST=EFLG)
	  IF (.NOT.EFLG) GOTO 5
        OPEN (UNIT=2,FILE=INFILE,STATUS='OLD')
c	READ(2,*) '.SPACE_GROUP_operators    r   12 (3f12.7)'
25	READ(2,'(A40)') string
	if (string(:1).ne.'.') goto 25
	do 50 j=1,4
C50	read(2,1151)(M1(I,J),I=1,3)
50	read(2,*)(M1(I,J),I=1,3)
1151	FORMAT(3f14.7)
	CLOSE(UNIT=2)
	MM=3
	NN=4
	CALL MATPRN (M1,MM,NN)

c 	read SECOND operator into m2
6	WRITE(6,*) 'ENTER NAME OF left hand (rot2,t2) rotation-trans op (.o):'
	READ(5,1200) INFILE
	IF (INFILE(:5).EQ.'     ') STOP
	  INQUIRE(FILE=INFILE, NAME=FILENAME,EXIST=EFLG)
	  IF (.NOT.EFLG) GOTO 6
7       OPEN (UNIT=2,FILE=INFILE,STATUS='OLD')
55	READ(2,'(A40)') string
	if (string(:1).ne.'.') goto 55
	do 60 j=1,4
60	read(2,*)(M2(I,J),I=1,3)
	CLOSE(UNIT=2)
	write(6,*)
	CALL MATPRN (M2,MM,NN)

c	Multiply m1 matrix and transl by m2 matrix (first 3 cols), put im M3

	do 100 i=1,3
	do 100 j=1,4
	m3(i,j)=0
	do 100 k=1,3
	m3(i,j)=m3(i,j)+m2(i,k)*m1(k,j) 
C	WRITE(6,*)I,J,K,M2(I,K),M1(K,J),M3(I,J)
100	CONTINUE
c	now add the translation of the second operator:
	do 110 i=1,3
110	m3(i,4)=m3(i,4)+m2(i,4)

	WRITE(6,*)' RESULTING NCS OPERATOR:'
	WRITE(6,*)' ROTATE MATRIX'
	DO 1100 I=1,3
1100	WRITE(6,9) (M3(I,J),J=1,3)
	WRITE(6,*)' TRANSLATE '
    	WRITE(6,9) (M3(I,4),I=1,3)
9	FORMAT(' ',3F14.6)



8	WRITE(6,*) 'ENTER NAME OF left hand (rot2,t2) rotation-trans op (.o):'
	READ(5,1200) INFILE
	IF (INFILE(:5).EQ.'     ') goto 1120
	  INQUIRE(FILE=INFILE, NAME=FILENAME,EXIST=EFLG)
	  IF (.NOT.EFLG) GOTO 8
c	copy m3 to m1
	do 102 i=1,3
	do 102 j=1,4
102	m1(i,j)=m3(i,j)
	goto 7

c	calculate product of diagonal elements (similarity to identity)
1120	x=1.
	do 1122 i=1,3
1122	x=x*m3(i,i)
	write(6,*)' diagonal product:',x

	write(6,*)' this operator will be saved as an O datablock RT.O'
	WRITE(6,*)'  Note this operator is for orthogonal coord, A translation'
c	WRITE(6,*)'  You may have to convert for your cell.    '
	WRITE(6,*) 'ENTER NAME FOR output of resulting rot-trans op (.o):'
	READ(5,1200) filename
	open (unit=3,name=FILENAME,status='unknown')
	WRITE(3,1149) 
1149	format('.lsq_rt_prod    r 12 (3f14.7) ')	
	do 1150 j=1,4
1150	write(3,1151)(M3(I,J),I=1,3)
	END

	SUBROUTINE MATPRN(ARRAY,M,N)
	REAL*8 ARRAY(3,4)
	WRITE(6,*)'M, N=',M,N
	DO 10 I=1,M
10	WRITE (6,50) (ARRAY(I,J),J=1,N)
50	format(8f10.5)
  	RETURN
        END