PRO SPHE_ROT_ANGLE,rok,mesic,den,tfirst,P0,drot,coel_stat,coel=COEL
;
; Calculation of the rotation angle DROT based on the daily motion of the Sun and
; on the position angle of the solar axis P0.
; INPUTS: rok = year, mesic = month, den = day, tfirst = UT of scan start, P0
; KEYWORD: coel - if set, the coelostat position is set according to the COEL value
; OUTPUT: drot = angle of rotation of the scanned spectroheliogram (used in DISKCOR.PRO)
;    to put the meridian in vertical direction (degrees)
;         coel_stat - info about coelostat position (string)
;
; Copyright M. Klvana and M. Sobotka, ASU AVCR Ondrejov, Czech Rep., 2008


srok=strtrim(rok,2)
srok=strmid(srok,2,2)
if mesic LT 10 then smes='0'+strtrim(mesic,2) else smes=strtrim(mesic,2)
if den LT 10 then sden='0'+strtrim(den,2) else sden=strtrim(den,2)

q=smes+sden
datum=srok+smes+sden
hodina=fix(tfirst/3600.)    ; UT in hours

if  ((q GE '1213') and (q  LT  '1232'))  then DCH=8.0
if  ((q GE '0101') and (q  LT  '0108')) or  ((q GE '1206') and (q  LT '1213')) then DCH=8.1
if  ((q GE '0108') and (q  LT  '0114')) or  ((q GE '1130') and (q  LT '1206')) then DCH=8.2
if  ((q GE '0114') and (q  LT  '0118')) or  ((q GE '1126') and (q  LT '1130')) then DCH=8.3
if  ((q GE '0118') and (q  LT  '0122')) or  ((q GE '1122') and (q  LT '1126')) then DCH=8.4
if  ((q GE '0122') and (q  LT  '0126')) or  ((q GE '1118') and (q  LT '1122')) then DCH=8.5
if  ((q GE '0126') and (q  LT  '0130')) or  ((q GE '1114') and (q  LT '1118')) then DCH=8.6
if  ((q GE '0130') and (q  LT  '0202')) or  ((q GE '1110') and (q  LT '1114')) then DCH=8.7
if  ((q GE '0202') and (q  LT  '0206')) or  ((q GE '1107') and (q  LT '1110')) then DCH=8.8
if  ((q GE '0206') and (q  LT  '0210')) or  ((q GE '1103') and (q  LT '1107')) then DCH=8.9
if  ((q GE '0210') and (q  LT  '0214')) or  ((q GE '1030') and (q  LT '1103')) then DCH=9.0
if  ((q GE '0214') and (q  LT  '0218')) or  ((q GE '1025') and (q  LT '1030')) then DCH=9.1
if  ((q GE '0218') and (q  LT  '0223')) or  ((q GE '1020') and (q  LT '1025')) then DCH=9.2
if  ((q GE '0223') and (q  LT  '0301')) or  ((q GE '1015') and (q  LT '1020')) then DCH=9.3
if  ((q GE '0301') and (q  LT  '0308')) or  ((q GE '1007') and (q  LT '1015')) then DCH=9.4
if  ((q GE '0308') and (q  LT  '1007')) then DCH=9.5


if hodina LE 13 then coelpos='E' else coelpos='W'

;in summer period is the coelostat in the central position C between E and W

;for new measurements: (summer period is from 20.03. to 20.09.)
if  ((q GE '0320') and (q LT '0920')) then coelpos='C'

; for old measurements:
if  ((datum GE '040319') and (datum LT '040921')) then coelpos='C'
if  ((datum GE '050315') and (datum LT '050918')) then coelpos='C'
if  ((datum GE '060322') and (datum LT '060920')) then coelpos='C'

;possible correction
if keyword_set(coel) then if coel eq 'E' or coel eq 'C' or coel eq 'W' then coelpos=coel
print,'Coelostat position:  ',coelpos
coel_stat=coelpos

if coelpos EQ 'C' then DCHR=0
if coelpos EQ 'E' then DCHR=(-DCH)
if coelpos EQ 'W' then DCHR=DCH

drot = DCHR+P0          ;P0 - position angle of the solar axis
         ;DROT in degrees
         ;P0 - uhel mezi osou rotace Slunce a kolmici k dennimu chodu
                        ;mereny od N konce kolmice smerem k ose rotace,
         ;smerem na je E kladny

print,'Rotation angle = ',drot,'  degrees'

END