PROGRAM SOLARM
C	FINDS COMPONENTS OF SOLAR MOTION USING INPUTTED COORDINATES
C       AND VELOCITIES OF A GROUP OF STARS
C	ENTER DATA IN UNITS OF RA, DEC, & RADIAL VELOCITY
C	RA AND DEC IN DECIMAL UNITS OF DEGREES
C	TAKES INPUT FILE LABELLED solarm.dat
	REAL  RA,DEC,VR
	REAL  A1,A2,A3,B1,B2,B3,C1,C2,C3,D1,D2,D3,E1,E2,E3,Q
	REAL  DET,DETX,DETY,DETZ,X,Y,Z,DETDX,DETDY,DETDZ,DX,DY,DZ
	REAL  DETEX,DETEY,DETEZ,EX,EY,EZ,F1,F2,F3
	DATA  A1,A2,A3,B1,B2,B3 / 0.0,0.0,0.0,0.0,0.0,0.0 /
	DATA  C1,C2,C3,D1,D2,D3 / 0.0,0.0,0.0,0.0,0.0,0.0 /
	DATA  E1,E2,E3,I / 0.0,0.0,0.0,0 /
	DIMENSION  A(100),D(100),V(100),P(100)
	OPEN (UNIT=7, FILE='solarm.dat', STATUS='OLD')
10	READ (7,*,END=12) RA, DEC, VR
	I = I + 1
	A(I) = RA
	D(I) = DEC
	V(I) = VR
	A1 = A1 + COSD(A(I))*COSD(A(I))*COSD(D(I))*COSD(D(I))
	A2 = A2 + COSD(A(I))*SIND(A(I))*COSD(D(I))*COSD(D(I))
	A3 = A3 + COSD(A(I))*COSD(D(I))*SIND(D(I))
	B1 = A2
	B2 = B2 + SIND(A(I))*SIND(A(I))*COSD(D(I))*COSD(D(I))
	B3 = B3 + SIND(A(I))*COSD(D(I))*SIND(D(I))
	C1 = A3
	C2 = B3
	C3 = C3 + SIND(D(I))*SIND(D(I))
	D1 = D1 - V(I)*COSD(A(I))*COSD(D(I))
	D2 = D2 - V(I)*SIND(A(I))*COSD(D(I))
	D3 = D3 - V(I)*SIND(D(I)) 
	GOTO 10
12	DET = A1*(B2*C3 - B3*C2) - A2*(B1*C3 - B3*C1) + A3*(B1*C2 - B2*C1)
	DETX = D1*(B2*C3 - B3*C2) - A2*(D2*C3 - B3*D3) + A3*(D2*C2 - B2*D3)
	DETY = A1*(D2*C3 - B3*D3) - D1*(B1*C3 - B3*C1) + A3*(B1*D3 - D2*C1)
	DETZ = A1*(B2*D3 - D2*C2) - A2*(B1*D3 - D2*C1) + D1*(B1*C2 - B2*C1)
	X = DETX/DET
	Y = DETY/DET
	Z = DETZ/DET
	DO 13 J = 1,I
	Q = Z*SIND(D(J))
	P(J) = 0.0 - X*COSD(A(J))*COSD(D(J)) - Y*SIND(A(J))*COSD(D(J)) - Q
	E1 = E1 - ABS(V(J) - P(J))*COSD(A(J))*COSD(D(J))
	E2 = E2 - ABS(V(J) - P(J))*SIND(A(J))*COSD(D(J))
	E3 = E3 - ABS(V(J) - P(J))*SIND(D(J))
13	CONTINUE
	F1 = E1/REAL(I - 1)
	F2 = E2/REAL(I - 1)
	F3 = E3/REAL(I - 1)
14	DETDX = E1*(B2*C3 - B3*C2) - A2*(E2*C3 - B3*E3) + A3*(E2*C2 - B2*E3)
	DETDY = A1*(E2*C3 - B3*E3) - E1*(B1*C3 - B3*C1) + A3*(B1*E3 - E2*C1)
	DETDZ = A1*(B2*E3 - E2*C2) - A2*(B1*E3 - E2*C1) + E1*(B1*C2 - B2*C1)
	DETEX = F1*(B2*C3 - B3*C2) - A2*(F2*C3 - B3*F3) + A3*(F2*C2 - B2*F3)
	DETEY = A1*(F2*C3 - B3*F3) - F2*(B1*C3 - B3*C1) + A3*(B1*F3 - F2*C1)
	DETEZ = A1*(B2*F3 - F2*C2) - A2*(B1*F3 - F2*C1) + F1*(B1*C2 - B2*C1)
	DX = DETDX/DET
	DY = DETDY/DET
	DZ = DETDZ/DET
	EX = DETEX/DET
	EY = DETEY/DET
	EZ = DETEZ/DET
	WRITE (6,15) X, DX, EX
15	FORMAT (1X,'X = ',F13.5,' DX = +-',F13.5,' SE = +-',F13.5)
	WRITE (6,16) Y, DY, EY
16	FORMAT (1X,'Y = ',F13.5,' DY = +-',F13.5,' SE = +-',F13.5)
	WRITE (6,17) Z, DZ, EZ
17	FORMAT (1X,'Z = ',F13.5,' DZ = +-',F13.5,' SE = +-',F13.5)
	STOP
	END