c This function returns the julian date of easter at 12 noon for c year in the gregorian calendar. c (see also: http://aa.usno.navy.mil/AA/) real*8 function easter (year) integer year ia = mod (year, 19) ib = year / 100 ic = year - ib * 100 id = ib / 4 ie = ib - id * 4 if = (ib + 8) / 25 ig = (ib - if + 1) / 3 ih = mod (ia * 19 + ib - id - ig + 15 , 30) ii = ic / 4 ik = ic - ii * 4 il = mod (32 + 2 * ie + 2 * ii - ih - ik , 7) im = (ia + 11 * ih + 22 * il) / 451 iq = ih + il - 7 * im + 22 number_of_leap_years = year / 4 + year / 400 - year / 100 easter = 1721119. d0 + year * 365.d0 . + number_of_leap_years + iq return end