################## Defining new types ################## Defining new types is done by implementing new reading methods in the corresponding modules (``pygeodyn/inout/priors.py`` for priors and ``pygeodyn/inout/observations.py`` for observations). .. warning:: For the modifications in the source code to be taken into account, the package must have been installed in editable mode (:code:`pip install -e`) ! See the `installation section`_ for more info. Priors ====== For examples, you can refer to the priors embarked with the program (that are outputs from geodynamo simulations): * Coupled-Earth [AFF13]_: ``0path`` * Midpath [AGF17]_: ``50path`` * 71% path [AG21]_: ``71path`` * 100% path [A23]_: ``100path`` The reading methods are found in `pygeodyn/inout/priors.py `_ . To use you own data, you must implement your reading method in the same file and design it with the same format. You can mix up priors (like ``71path`` with ``100path``). An example of the mix of 71path (long time series) with 100path (short time series) can be found in read_100_path (`pygeodyn/inout/priors.py `_). Namely: .. code-block:: python # The method name must have this name format and signature def read_(data_directory, prior_type, dt_sampling, measures=None): # Create containers to store priors MF, U, ER, times, dt_samp, tag = [], [], [], [], [], [] # Read every prior... # If mix of prior then THE LONGEST TIME SERIES PRIOR MUST BE READ AND STORED FIRST # if no time in prior then times = None and dt_samp = None # ...append prior to containers # ...and returns the containers return MF, U, ER, times, dt_samp, tag Now, all is left is to set ``prior_type`` in the config file to the name ```` that you chose for it (and ``prior_data_dir`` to the directory where the prior data are stored). Then, try to launch a computation and check if any errors are raised (``NotImplementedError`` if your prior is not associated to a reading function, ``IOError`` if the files were not found...). If not, your prior data will be used as average data in the AR-1 process and to compute the covariance matrices. Observations ============ Examples of valid reading methods are given for the embarked observations: * SOLA [HF19]_: ``sola_hdf5`` * CHAOS [FKO20]_: ``chaos_hdf5`` * COVOBS [GJF13]_, [GBF15]_: ``covobs_hdf5`` * KALMAG [JMS22]_: ``kalmag_hdf5`` * Ground and Virtual observatories [BHF18]_, [Hammer18]_: ``go_vo`` The strategy to use another type of observations follows the same paradigm as for priors: you must define a new observation type by implementing its reading method in `pygeodyn/inout/observations.py `_. In more details, the format must be: .. code-block:: python # The method name must have this name format and signature def build__observations(cfg, nb_realisations, measure_type): # Read your observation data for measure_type='MF' or 'SV'... # ...and returns a dict of Observation objects with numpy.datetime64 objects as keys return observations The implementation is trickier in this case and several points need to be highlighted: * Your new reading method must take the measure_type as arg that will be equal to either 'MF' (for magnetic field observations) or 'SV' (for secular variation observations). * The returned object must be a dictionary with `numpy.datetime64 `_ objects as keys that correspond to the date of the observations * The said dictionary must contain `Observation objects `_ that are a handy way to store all the observation parameters: * The data :math:`Y` as 2D numpy.ndarrays (nb_realisations x nb_observed_coefs) * The observation operator :math:`H` (nb_observed_coefs x nb_coefs) that verifies :math:`Y = HX` where :math:`X` is a vector containing the spectral coefficients of the magnetic field or secular variation (size: nb_coefs). * The observation errors matrices :math:`R` (nb_observed_coefs x nb_observed_coefs) For reference, here is the basic way to create an *Observation* object and the *observations* dictionary: .. code-block:: python observations = {} # ... # Assuming you have extracted all the aforementioned parameters observation_at_t = Observation(obs_data, H, R) #Let say you consider analyses at time 1997, 1997.5, 1998, 1998.5 #Then to add an observation, you must put the analysis index as the observation attribute #for 1997 observations["0"] = observation_at_t #for 1997.5 observations["1"] = observation_at_t # ... Do this for all observation dates and return the observations dict After having implemented your new method, set ``obs_type`` in the config file to the name ```` that you chose for it (and ``obs_dir`` to the directory where the data are stored). Then, try to launch a computation and check if any errors are raised or logged (if no reading function are found, you will get a warning in the logs saying that all analyses will be skipped). If not, your observation data will be used as observations in the analysis step of the algorithm. Example ======= .. admonition:: Work in progress Write an example with all the steps References ========== .. [AFF13] Aubert, J., Finlay, C. C. & Fournier, A. Bottom-up control of geomagnetic secular variation by the Earth’s inner core, *Nature*, 502, 219–223, 2013. doi:10.1038/nature12574 .. [AGF17] Aubert, J., Gastine, T. & Fournier, A. Spherical convective dynamos in the rapidly rotating asymptotic regime, *Journal of Fluid Mechanics*, 813, 558–593, 2017. doi:10.1017/jfm.2016.789 .. [AG21] Aubert, J. & Gillet, N. The interplay of fast waves and slow convection in geodynamo simulations nearing Earth’s core conditions, *Geophysical Journal International*, 225, 1854–1873, 2021. doi:10.1093/gji/ggab054 .. [A23] Aubert, J.: State and evolution of the geodynamo from numerical models reaching the physical conditions of Earth’s core, Geophys. J. Int. 235, 468-487, 2023, doi: 10.1093/gji/ggad229 .. [BHF18] Barrois, O., Hammer, M. D., Finlay, C. C., Martin, Y. & Gillet, N. Assimilation of ground and satellite magnetic measurements: inference of core surface magnetic and velocity field changes, *Geophysical Journal International*, 2018. doi:10.1093/gji/ggy297 .. [GBF15] Gillet, N., Barrois, O. & Finlay, C. C. Stochastic forecasting of the geomagnetic field from the COV-OBS.x1 geomagnetic field model, and candidate models for IGRF-12, *Earth, Planets and Space*, 67, 2015. doi:10.1186/s40623-015-0225-z .. [GJF13] Gillet, N., Jault, D., Finlay, C. C. & Olsen, N. Stochastic modeling of the Earth’s magnetic field: Inversion for covariances over the observatory era, *Geochemistry, Geophysics, Geosystems*, 14, 766–786, 2013. doi:10.1002/ggge.20041 .. [Hammer18] Hammer, M. D. Local Estimation of the Earth’s Core Magnetic Field. (Technical University of Denmark (DTU), 2018). .. [IGF23] Istas, M., Gillet, N., Finlay, C. C., Hammer, M. D., & Huder, L. Transient core surface dynamics from ground and satellite geomagnetic data, *Geophysical Journal International*, 233, 1890–1915, 2023. doi:10.1093/gji/ggad039 .. [FKO20] Finlay CC, Kloss C, Olsen N, Hammer MD, Tøffner-Clausen L, Grayver A, Kuvshinov A. The CHAOS-7 geomagnetic field model and observed changes in the South Atlantic Anomaly. Earth Planets Space. 2020;72(1):156. doi: 10.1186/s40623-020-01252-9 .. [JMS22] Julien, B., Matthias, H., Saynisch-Wagner, J. et al. Kalmag: a high spatio-temporal model of the geomagnetic field. Earth Planets Space 74, 139 (2022). https://doi.org/10.1186/s40623-022-01692-5 .. [HF19] Magnus D Hammer, Christopher C Finlay, Local averages of the core–mantle boundary magnetic field from satellite observations, Geophysical Journal International, Volume 216, Issue 3, March 2019, Pages 1901–1918, https://doi.org/10.1093/gji/ggy515 .. _installation section: https://geodynamo.gricad-pages.univ-grenoble-alpes.fr/pygeodyn/install.html#install-the-pygeodyn-package-from-sources