vallis / libstempo

libstempo — a Python wrapper for tempo2
MIT License
29 stars 34 forks source link

Allow time array to be passed directly to the tempopulsar class #37

Closed mattpitkin closed 3 years ago

mattpitkin commented 3 years ago

Inspired by #35, I thought it would be useful to be able to pass arbitrary arrays of TOAs directly to the tempopulsar object rather than through the fakepulsar function. The implementation here means that a temporary .tim file does not have to be written to file if you want to create a pulsar object.

This PR also exposes the phase values and provides a method to convert residuals to phase residuals.

jellis18 commented 3 years ago

Could you maybe add an example of this to one of the example notebooks?

mattpitkin commented 3 years ago

Could you maybe add an example of this to one of the example notebooks?

Yep, will do.

vallis commented 3 years ago

This is a very extensive change that effectively replicates the .tim-file reading functionality of tempo2. Do we then need a plan to keep this in sync with tempo2? How are we sure that we're assigning/allocating everything that needs to be? Is this kind of functionality best left for PINT (assuming it supports making fake pulsars, which I don't know)?

More technically, do inputtoas, inputtoaerrs, inputobservatory, inputobsfreqs need to have public interfaces?

This same PR has also updates to residuals() (which could use some comments before I can review it) and the new phaseresiduals() and phase(), which could perhaps be implemented as one-liners (or almost) by calling self.residuals(**kwargs).

mattpitkin commented 3 years ago

More technically, do inputtoas, inputtoaerrs, inputobservatory, inputobsfreqs need to have public interfaces?

Good point. I've now switched these to set the values using private methods rather than properties.

mattpitkin commented 3 years ago

Is this kind of functionality best left for PINT (assuming it supports making fake pulsars, which I don't know)?

Just as a bit a background for the main reason I'm doing this: I want to be able to use TEMPO2 to calculate the pulsar phases for use in my code that searches for GWs from pulsars (which I'm currently rewriting in Python). This code currently uses custom solar system and binary barycentring routines within LALSuite (the SSB routines were written quite a while ago by Curt, based on those in TEMPO), but it has many limitations when trying to use pulsar par files fit with more "modern" parameters. PINT would be great to use, as it's already natively Python, but has the one major downside - it still can't handle TCB-based pulsar ephemerides which are by default produced by TEMPO2, and getting this ability into PINT is non-trivial. libstempo therefore seems like a very nice way to go. For my purpose, I need to calculate phases for a large number of time values, so writing these to a .tim file each time would seem slow and inefficient, hence the changes suggested rather than making use of fakepulsar.

This is a very extensive change that effectively replicates the .tim-file reading functionality of tempo2. Do we then need a plan to keep this in sync with tempo2? How are we sure that we're assigning/allocating everything that needs to be?

As above, I've essentially added this for my own reasons and don't know if other people would want to make use of it, so I'm happy to maintain this. I haven't added in the ability to add all the things that can be input in a .tim file, for instance, I don't allow the setting of any flags or EFAC/EQUAD values. But, from my testing I get identical results between fakepulsar and the direct input method. I'll add a comparison directly to the test suite, so that any issues can be caught if/when something changes in TEMPO2.

mattpitkin commented 3 years ago

I've added a test file for both the fakepulsar function and passing TOAs via the tempopulsar object.

mattpitkin commented 3 years ago

@jellis18 I've added an example to the libstempo-toasim-demo.ipynb notebook.

mattpitkin commented 3 years ago

This same PR has also updates to residuals() (which could use some comments before I can review it) and the new phaseresiduals() and phase(), which could perhaps be implemented as one-liners (or almost) by calling self.residuals(**kwargs).

I've switched phaseresiduals() and phase() to just take **kwargs as an argument and pass that on to residuals() and added some (hopefully) explanatory comments about the updates to residuals (see also the discussion here).

mattpitkin commented 3 years ago

@vallis @jellis18 I've fixed the CI tests that were previously failing. Is there anything else you'd like to see tested/documented with these changes?

vallis commented 3 years ago

Should I make a release now or wait for further development?

mattpitkin commented 3 years ago

Thanks very much @vallis. I'd appreciate a new release now, but can wait for other developments to be added if necessary.

mattpitkin commented 3 years ago

@vallis I don't see any other PR or issues that would suggest other updates being added soon, so would you be able to create a new release now? Just bumping to v2.4.2 would be fine.