search

SWxTREC / pymsis

Python interface to the NRLMSIS codes
MIT License
23 stars 5 forks source link

Daily F10.7 is off by one day in get_f107_ap #54

Closed rogervarney closed 2 months ago

rogervarney commented 2 months ago

It looks like the utility to get F10.7 is off by one day. This minimum working example

import numpy as np from pymsis import msis from pymsis import utils as msisutils

date = np.datetime64("2016-03-04T02:00") f107,f107a,ap = msisutils.get_f107_ap(date) print(f107)

Prints 98.7, which would the correct value for 2016-03-03. The correct value for 2016-03-04 should be 100.5. I've checked against both NGDC (https://www.ngdc.noaa.gov/stp/space-weather/solar-data/solar-features/solar-radio/noontime-flux/penticton/penticton_observed/listings/listing_drao_noontime-flux-observed_daily.txt) and Celestrak (https://celestrak.org/SpaceData/).

The Ap outputs look correct. It's just F10.7 that is off by one day.

greglucas commented 2 months ago

This is one area I was following documentation, but I wasn't completely sure why they defined it that way. I haven't read the MSIS papers in a while now. But, I defined this get_f107_ap() to be MSIS specific and return the values that MSIS needs for input, which is the previous day's value.

https://github.com/SWxTREC/pymsis/blob/d30ab2744c5e6033ad0ef5d197021571fa7ef1a6/pymsis/utils.py#L191-L192

Here is a link to the MSIS source comments for that in GitHub: https://github.com/jacobwilliams/NRLMSIS2.1/blob/912414461c84440f157a8761c1f2ca1efb03c522/msis_gtd8d.F90#L34-L47

Let me know if you think there is an error still, or if we should try to document this better somehow.

rogervarney commented 2 months ago

Wow you are right. Seems like I have been running MSIS incorrectly for years.

The appendix of Hedin 1987 in explaining equation A22 clearly says "F10.7 = 10.7 cm flux on previous day, 10^-22 W m^-2 Hz^-1, \bar{F}10.7 is the average F10.7 over three solar rotations (81 days) centered on the required day, ..."