C+
C NAME:
C	Time2PAnglePole
C PURPOSE:
C	Calculate the position angle of the solar poles
C CATEGORY:
C	Astronomy: celestial physics
C CALLING SEQUENCE:
	double precision function Time2PAnglePole(id,tt,lng)
C INPUTS:
C	ID	integer			= 1,-1 : Get PA relative to Equatorial North of Earth
C					= 2,-2 : Get PA relative to Ecliptic North
C					ID>0 gets the PA of the north pole
C					ID<0 gets the PA of the south pole
C	tt 	integer			time
C	lng	double precision	Ecliptic longitude of Sun as viewed from
C					observer (i.e. topocentric longitude)
C OUTPUTS:
C	Time2PAnglePole
C		double precision	Position angle in degrees
C INCLUDE:
	include		'sun.h'
C CALLS:
C	Time2YDoy, Time2Day8
C PROCEDURE:
C >	For given time tt and ecliptic longitude of the Sun (as seen
C	from an observer in the ecliptic plane) the position of angle
C	of the solar pole is calculated, either relative to ecliptic or to
C	equatiorial North.
C	See Spherical Astronomy pp. 430-433 by Green, Robin M., Cambridge UP (1985)
C	(The position angle is measured counterclockwise, i.e. toward the east,
C	from North)
C >	PA(south pole) = PA(north pole)+180
C MODIFICATION HISTORY:
C	JUL-1992, Tom Davidson, Paul Hick (UCSD/CASS)
C	JAN-2004, Paul Hick (UCSD/CASS)
C		Introduced 2-elements integer standard time
C-
	    integer		id
	    integer		tt(2)
	    double precision	lng

	double precision	E		/23.44d0/	! Obliquity of the ecliptic
	integer			yr
	integer			doy(2)

	double precision	T
	double precision	O

	double precision	dcosd
	double precision	dtand
	double precision	datand

	call Time2YDoy(0,tt,yr,doy)
	call Time2Day8(0,doy,T)

	T = yr-1850+(T-1.0d0)/365.25d0		! # years since 1850
	O = 73.0d0+40.0d0/60.0d0 + 0.5025d0/36.0d0*T ! Longitude of ascending node of solar 
						! equator on ecliptic
	Time2PAnglePole = -datand(dtand(SUN__I)*dcosd(0.0d0-lng))
	if (abs(id) .eq. 1) Time2PAnglePole = Time2PAnglePole-datand(dcosd(lng)*dtand(E))

	if (id .lt. 0) Time2PAnglePole = Time2PAnglePole+180.0d0 ! South pole

	return
	end