Recent Changes - Search:

edit SideBar

Development

Appendix A: Developing the software

The plug-in features of tempo2 implies that users can add to the functionality of \textsc{tempo2} by producing their own plug-in. In this section, we first discuss how such plug-ins can be produced. We also provide a few notes for users who are interested in developing the actual main source code and become part of the tempo2 development team.

A.1 Creating a plug-in

Plug-ins are relatively easy to create in C or C++. These plug-ins may use other libraries such as pgplot, tcl/tk or qt. As described above, plug-ins can be created to change the output format or as a graphical interface.

When creating plugins with command-line arguments, do notice that the standard tempo2 command line arguments are already in use and duplication may cause your plugin to behave unexpectedly.

A.1.1 A new output format

When an "output-format" plug-in is used, the tempo2 software reads in the specified parameter and observation files, carries out the standard techniques to obtain pre- and post-fit timing residuals and their corresponding parameter values and then passes control to the output-format software to display the results. An example output format (in C) is as follows:

 #include <stdio.h>
 #include <stdlib.h>
 #include <string.h>
 #include <math.h>
 #include "tempo2.h"  /* Essential for a tempo2 plugin */

 extern "C" int tempoOutput(pulsar *psr,int npsr,char timFile[][MAX_FILELEN],char parFile[][MAX_FILELEN]) 
 {  
   int i;
   printf("Pulsar name = %s\n",psr[0].name);
   printf("Pulsar spin-frequency = %Lg\n",psr[0].param[param_f].val[0]);
   for (i=0;i<psr[0].nobs;i++)
     printf("Residual = %.14f\n",(double)psr[0].obsn[i].residual);
 }

This plug-in displays the pulsar's name (obtained from the parameter file), its post-fit spin-frequency and the post-fit timing residuals on the screen. The software is built around a "structure" called "pulsar" which is defined in "tempo2.h". The most important definitions in this structure are

   char *name;                  /* The pulsar name */
   char *binaryModel;           /* The binary model (e.g. T2/BT/ELL1 ...) */
   double ne_sw;                /* The electron density at 1AU due to the solar wind */
   int    fitMode;              /* = 1 if fitting with errors, 0 = not fitting with errors */        
   int    nobs;                 /* Number of observations */
   int    nits;                 /* Number of iterations used for the fit */
   int    ipm;                  /* = 1 if the interplanetary DM correction is used, 0 otherwise */

An array of pulsars is defined and therefore, in order to display the name of the first three pulsars analysed:

   printf("pulsar 1 = %s\n",psr[0].name);
   printf("pulsar 2 = %s\n",psr[1].name);
   printf("pulsar 3 = %s\n",psr[2].name);

The structure also includes an array of parameters ("param") which are defined as

   param_raj     : Right ascension of the source
   param_decj    : Declination of the source
   param_f       : Pulse frequency (and derivatives)
   param_pepoch  : Epoch of period determination
   param_posepoch: Epoch of position determination
   param_dmepoch : Epoch of DM determination
   param_dm      : Dispersion measure
   param_pmra    : Proper motion in right ascension
   param_pmdec   : Proper motion in declination
   param_px      : Parallax
   param_sini    : Sine of the orbital inclination angle
   param_pb      : Orbital period
   param_t0      : Epoch of periastron passage
   param_a1      : Projected semi-major axis of the binary orbit
   param_om      : Angle of periastron
   param_pmrv    : Radial proper motion
   param_ecc     : Eccentricity of the binary orbit
   param_edot    : First time derivative of the eccentricity
   param_xpbdot  : Part of pbdot unexplained by GR (DDGR model)
   param_pbdot   : First time derivative of the orbital period
   param_a1dot   : First time derivative of the projected semi-major axis
   param_omdot   : First time derivative of the angle of periastron (periastron advance) 
   param_tasc    : Epoch of ascending node passage
   param_eps1    : Eccentricity parameter
   param_eps2    : Eccentricity parameter
   param_m2      : Companion mass
   param_gamma   : Gravitational redshift parameter
   param_mtot    : Total mass of the binary system
   param_glep    : 
   param_glph    : 
   param_glf0    : 
   param_glf1    : First time derivative of glf0
   param_glf0d   :
   param_gltd    : 
   param_start   : Start MJD of this data set
   param_finish  : Final MJD of this data set
   param_track   : 
   param_bp      : 
   param_bpp     : 
   param_tzrmjd  : 
   param_tzrfrq  : 
   param_fdd     : 
   param_dr      : 
   param_dtheta  :
   param_tspan   : Time span of the data set
   param_bpjep   : 
   param_bpjph   :
   param_bpja1   :
   param_bpjec   : 
   param_bpjom   : 
   param_bpjpb   :
   param_wave_om : Frequency of the fitted wave (for pre-whitening)
   param_kom     : Longitude of the ascending node of the binary orbit
   param_kin     : Inclination angle of the binary orbit
   param_shapmax : Maximum of the Shapiro delay (-ln(1-s))
   param_dth     : 
   param_a0      : 
   param_b0      : 
   param_xomdot  : Part of omdot unexplained by GR (DDGR)
   param_afac    : 
   param_eps1dot : First time derivative of eps1
   param_eps2dot : First time derivative of eps2
   param_tres    : Residual rms

For each parameter, the user can obtain

   char **label;              /* Label about this parameter                         */
   char **shortlabel;         /* Label about this parameter without units           */
   longdouble *val;           /* Value of parameter                                 */
   longdouble *err;           /* Uncertainty on parameter value (Notice the fitting routine only calculates 
                                 this parameter up to double precision.)
   int  *fitFlag;             /* = 1 if fitting required                            */
   int  *paramSet;            /* = 1 if parameter has been set                      */
   longdouble *prefit;        /* Pre-fit value of the parameter                     */
   longdouble *prefitErr;     /* Pre-fit value of the uncertainty                   */
   int aSize;                 /* Number of elements in the array for this parameter */

Each parameter (except those for which the derivatives were explicitly mentioned before) is stored as an array with each element of the array typically storing time derivatives of the parameter. For instance,to obtain the value of the spin-frequency and its first two derivatives:

   printf("values = Lg %Lg\n",psr[0].param[param_f].val[0],
             psr[0].param[param_f].val[1],psr[0].param[param_f].val[2]);

(Note: 'C' requires the symbol 'Lf' to display values in 'long double' precision.)

The "pulsar" structure also contains the set of observations and corresponding parameters. The observation structure contains (other less-common parameters are defined in the tempo2.h file)

  longdouble sat;                 /* Site arrival time (MJD)                                    
  longdouble bat;                 /* Barycentric arrival time (MJD)                             
  int deleted;                    /* = 1 if observation has been deleted                        
  longdouble prefitResidual;      /* Pre-fit residual                                                   
  longdouble residual;            /* residual (in sec)                                                  
  double      freq;               /* Frequency of observation (in MHz)                          
  double      freqSSB;            /* Frequency of observation in barycentric frame (in Hz)      
  double      toaErr;             /* Error on TOA (in us)                                       
  char        fname[MAX_FILELEN]; /* Name of data file giving TOA                               
  char        telID[100];         /* Telescope ID                                               
  longdouble correctionTT_TB;     /* Correction to TDB/TCB                                      
  double einsteinRate;            /* Derivative of correctionTT_TB                              
  longdouble correctionTT_Teph;   /* Correction to Teph                                         
  longdouble correctionUT1;       /* Correction from site TOA to UT1                            
  double shapiroDelaySun;         /* Shapiro Delay due to the Sun                               
  double shapiroDelayJupiter;     /* Shapiro Delay due to Jupiter                               
  double shapiroDelaySaturn;      /* Shapiro Delay due to Saturn                                
  double shapiroDelayUranus;      /* Shapiro Delay due to Uranus                                
  double shapiroDelayNeptune;     /* Shapiro Delay due to Neptune                               
  double troposphericDelay;       /* Delay due to neutral refraction in atmosphere              
  double tdis1;                   /* Interstellar dispersion measure delay                      
  double tdis2;                   /* Dispersion measure delay due to solar system               
  longdouble roemer;              /* Roemer delay                                               

For instance, in order to display the site arrival-time, the barycentric arrival-time, and the Solar Shapiro delay for each observation:

  for (i=0;i<psr[0].nobs;i++)
   printf("values = Lf %f\n",psr[0].obsn[i].sat,psr[0].obsn[i].bat,psr[0].obsn[i].shapiroDelaySun);

A.1.1 A new graphical interface

When a graphical interface is used, tempo2 checks the command-line arguments, initialises the memory and then passes control to the graphical interface. The reading of the parameter and observation files, calculating barycentric arrival times, obtaining residuals, fitting and displaying the output must all be done within the graphical interface. The following tempo2 functions are commonly called:

   readParfile(...)           /* Read the parameter file             */
   readTimfile(...)           /* Read the observation file           */
   preProcess(...)            /* Needs to be called after reading the
                                 parameter and observation files     */
   formBatsAll(...)           /* Forms the barycentric arrival times */
   formResiduals(...)         /* Forms the timing residuals          */
   doFit(...)                 /* Calls the fitting algorithms        */
   textOutput(...)            /* The standard output display         */

An example of a simple interface would be

 #include <stdio.h>
 #include <string.h>
 #include <stdlib.h>
 #include <math.h>
 #include "tempo2.h"

 using namespace std;

 /* The main function called from the TEMPO2 package is 'graphicalInterface' */
 /* Therefore this function is required in all plugins                       */
 extern "C" int graphicalInterface(int argc,char *argv[],pulsar *psr,int *npsr) 
 {
  char parFile[MAX_PSR][MAX_FILELEN];
  char timFile[MAX_PSR][MAX_FILELEN];
  int i;
  double globalParameter;

  *npsr = 1;  /* For a graphical interface that only shows results for one pulsar */

  printf("Graphical Interface: name\n");
  printf("Author:              author\n");
  printf("Version:             version number\n");

  /* Obtain the .par and the .tim file from the command line */
  if (argc==4) /* Only provided .tim name */
    {
      strcpy(timFile[0],argv[3]);
      strcpy(parFile[0],argv[3]);
      parFile[0][strlen(parFile[0])-3] = '\0';
      strcat(parFile[0],"par");
    }

  /* Obtain all parameters from the command line */
  for (i=2;i<argc;i++)
    {
      if (strcmp(argv[i],"-f")==0)
	{
	  strcpy(parFile[0],argv[i+1]); 
	  strcpy(timFile[0],argv[i+2]);
	}
    }

  readParfile(psr,parFile,timFile,*npsr); /* Load the parameters       */
  readTimfile(psr,timFile,*npsr); /* Load the arrival times    */
  preProcess(psr,*npsr,argc,argv);

  for (i=0;i<2;i++)                   /* Do two iterations for pre- and post-fit residuals*/
    {
      formBatsAll(psr,*npsr);         /* Form the barycentric arrival times */
      formResiduals(psr,*npsr,0.0);    /* Form the residuals                 */
      if (i==0) doFit(psr,*npsr,&globalParameter,0,0);   /* Do the fitting     */
      else textOutput(psr,*npsr,globalParameter,0,0,0,"");  /* Display the output */
    }

    /* Now you have the parameters and residuals */
    /* which can be displayed or plotted         */


  return 0;
 }

A.1.2: The main source code

The tempo2 software can be obtained from the ATNF CVS repository soft_atnf/tempo2 (please contact members of the ATNF pulsar group for more information). To install type

   make install

The plug-ins provided with the source code are stored in soft_atnf/tempo2/plugin.


Edit - History - Print - Recent Changes - Search
Page last modified on June 03, 2008, at 02:45 PM