معرفی شرکت ها
GalDynPsr-1.0.0
تبلیغات ما
مشتریان به طور فزاینده ای آنلاین هستند. تبلیغات می تواند به آنها کمک کند تا کسب و کار شما را پیدا کنند.
مشاهده بیشتر
تبلیغات ما
مشتریان به طور فزاینده ای آنلاین هستند. تبلیغات می تواند به آنها کمک کند تا کسب و کار شما را پیدا کنند.
مشاهده بیشتر
تبلیغات ما
مشتریان به طور فزاینده ای آنلاین هستند. تبلیغات می تواند به آنها کمک کند تا کسب و کار شما را پیدا کنند.
مشاهده بیشتر
تبلیغات ما
مشتریان به طور فزاینده ای آنلاین هستند. تبلیغات می تواند به آنها کمک کند تا کسب و کار شما را پیدا کنند.
مشاهده بیشتر
تبلیغات ما
مشتریان به طور فزاینده ای آنلاین هستند. تبلیغات می تواند به آنها کمک کند تا کسب و کار شما را پیدا کنند.
مشاهده بیشترتوضیحات
A practise package
| ویژگی | مقدار |
|---|---|
| سیستم عامل | OS Independent |
| نام فایل | GalDynPsr-1.0.0 |
| نام | GalDynPsr |
| نسخه کتابخانه | 1.0.0 |
| نگهدارنده | [] |
| ایمیل نگهدارنده | [] |
| نویسنده | Dhruv Pathak |
| ایمیل نویسنده | pathakdhruv9786@gmail.com |
| آدرس صفحه اصلی | https://github.com/pathakdhruv/GalDynPsr |
| آدرس اینترنتی | https://pypi.org/project/GalDynPsr/ |
| مجوز | New BSD |
# GalDynPsr Package
GalDynPsr is a package for calculating dynamical contributions in the period derivatives of pulsars in the Galactic field. It can calculate the fractional contributions or the excess terms, e.g. \dot{P}/P|_excess where P is either the orbital period or the spin period. Various dynamical contributions, including the Shklovskii effect (the contribution due to the proper motion) or due to the acceleration of the pulsar caused by the gravitational potential of the Galaxy, can be calculated. Various methods or models to perform these tasks are available in the package. Using the measured values of the periods and period derivatives, GalDynPsr can even compute the "intrinsic" values of the period derivatives, provided no other extra contribution exist.
Details on various dynamical effects and various models to estimate those are available in Pathak and Bagchi (arXiv:1712.06590).
Brief outline of usage of GalDynPsr is given below.
Questions and queries can be sent to Dhruv Pathak pathakdhruv9786 at gmail dot com
# 1) After installing numpy, scipy, and galpy successfully, install GalDynPsr
pip3 install GalDynPsr
### If wished, one can change the values of Rs (galactocentric cylindrical radius of the sun), Vs (rotational speed of the sun) etc in parameters.in file that can be found inside the GalDynPsr (installed directory).
## But remember that galpy also has these values defined in the file '$home/.galpyrc'. One can in principle change these values.
# However, the Milky Way potential in galpy was fitted with Rs = 8.0 kpc and vs = 220.0 in galpy
# 2) Import GalDynPsr in your code as
import GalDynPsr
# 3) In all of the following modules, the ordering of the arguments of the functions are very important.
# 4) Calculate excess terms of modelX ** Very important, model names are case sensitive!
GalDynPsr.modelX.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) and GalDynPsr.modelX.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd)
The function Expl() calculates the excess term due to the relative acceleration of the pulsar parallel to the Galactic plane and the function Exz() calculates the excess term due to the relative acceleration of the pulsar perpendicular to the Galactic plane. Total dynamical contribution will be the sum of above two terms. Both of the above functions takes the arguments as the Galactic longitude in degrees (say ldeg), uncertainty in the Galactic longitude in degrees (say sigl), the Galactic latitude in degrees (say bdeg), uncertainty in the Galactic latitude in degrees (say sigb), the distance of the pulsar from the solar system barycenter in kpc (say dkpc), and the uncertainty in the distance in kpc (say sigd); strictly in this order. One needs to assign the values of ldeg, sigl, bdeg, sigb, dkpc, and sigd before calling the above functions.
The models that do not involve galpy can return uncertainities in the excess terms following standard error propagation method, like
modelX.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) and modelX.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd)
where the function Errpl() returns the uncertainty in the excess terms due to the relative acceleration of the pulsar parallel to the Galactic plane and Errz() returns the uncertainty in the excess terms due to the relative acceleration of the pulsar perpendicular to the Galactic plane. The arguments are the same as before, and needs to be assigend values before using the functions.
Remeber to replace modelX by the model you want to use, like modelA, modelLa, etc. Model names are case sensitive!
# 5) Calculate the Shklovskii term
The Shklovskii term can be calculated as GalDynPsr.Shk.Exshk(dkpc, sigd, mu_alpha, sigmua, mu_delta, sigmud) and the error in the Shklovskii term as GalDynPsr.Shk.errShk(dkpc, sigd, mu_alpha, sigmua, mu_delta, sigmud)
where mu_alpha is the proper motion in the Right Ascesion direction (mas/yr), mu_delta is the proper motion in the Declination direction (mas/yr), sigmua is the uncertainty in mu_alpha, and sigmud is the uncertainty in mu_delta. dkpc is as usual the distance of the pulsar from the solar system barycenter in kpc and sigd is the uncertainty in dkpc. One needs to assign the values of these parameters first.
# 6) Print the basic parameters of the pulsars
GalDynPsr.read_parameters.Rskpc returns the Galactocentric cylindrical radius of the sun.
GalDynPsr.read_parameters.Vs returns the rotational speed of the sun.
GalDynPsr.read_parameters.Rpkpc(ldeg, sigl, bdeg, sigb, dkpc, sigd) returns the value of Galactocentric cylindrical radius of the pulsar in kpc and GalDynPsr.read_parameters.z(ldeg, sigl, bdeg, sigb, dkpc, sigd) returns the verical height of the pulsar from the Galactic plane. The meaning of the arguments are as usual:
One can get the uncertainty in the Galactocentric cylindrical radius of the pulsar in kpc using GalDynPsr.read_parameters.ErrRpkpc(ldeg, sigl, bdeg, sigb, dkpc, sigd) and the uncertainty in the verical height as GalDynPsr.read_parameters.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd).
Obviously, one first needs to assign the values of ldeg, sigl, bdeg, sigb, dkpc, and sigd.
# 7) For a pulsar in a globular cluster, one can read the cluster parameters, e.g. the Galactic longitude, latitude and the distance from the Harric catalogue (GC_harris.txt coming with GalDynPsr), and calculate crude values of the dynamical terms
One first needs to import the module to read Harris cataloge as:
import GalDynPsr.ReadGC as readgc
Then one can input the name of the cluster, e.g. readgc.gcinput("47Tuc")
Then the excess terms using modelX can be calculated as
sigl=0.000001
sigb=0.000001
sigd=0.000001
GalDynPsr.modelX.Expl(readgc.ldeg,sigl,readgc.bdeg,sigb,readgc.dkpc,sigd) # the parallel component
GalDynPsr.modelX.Exz(readgc.ldeg,sigl,readgc.bdeg,sigb,readgc.dkpc,sig) # the perpendicular component
Note that the Harris catalogue does not provide uncertainties in ldeg, bdeg and dkpc. That is why we need to define them before calculating the excess terms. As the pulsar parameters differ from the cluster parameters as in Harris catalogue, it is meaningless to calculate uncetrainities. So, one SHOULD simply call the above functions as:
GalDynPsr.modelX.Expl(readgc.ldeg,0.0,readgc.bdeg,0.0,readgc.dkpc,0.0) # the parallel component
GalDynPsr.modelX.Exz(readgc.ldeg,0.0,readgc.bdeg,0.0,readgc.dkpc,0.0) # the perpendicular component
# 8) models involving galpy does not return error. However, that can be done using simulation. See Example-2 below.
# 9) Calculation of instrinsic period. Ordering of the arguments is important.
First assign values of ldeg, sigl, bdeg, sigb, dkpc, sigd (meaning of the parameters are as usual) and calculate the excess terms as usual:
Ex_pl = GalDynPsr.modelA.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) # excess term parallel to the Galactic plane
Ex_z = GalDynPsr.modelA.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) # excess term perpendicular to the Galactic plane
errpl = GalDynPsr.modelA.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) # uncertainty in the parallel term
errz = GalDynPsr.modelA.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) # uncertainty in the parallel term
Ex_shk = GalDynPsr.Shk.Exshk(dkpc, sigd, mu_alpha, sigmua, mu_delta, sigmud) # excess term due to the Shklovskii effect
errshk = GalDynPsr.Shk.errShk(dkpc, sigd, mu_alpha, sigmua, mu_delta, sigmud) # uncertainties in the Shklovskii term
Now assign the values of the period 'P' in seconds, the uncertainy in period 'sigP' in seconds, the measured value of the period derivative 'Pdot_obs' in seconds/seconds and the uncertainty in the measured value of the period derivative 'sigPdot' in seconds/seconds; and calculate the dynamically caused values of the period derivatives (in seconds/seconds) as:
GalDynPsr.pdotint.PdotExpl(Ex_pl,P) # due to the parallel component of the acceleration
GalDynPsr.pdotint.PdotExz(Ex_z,P) # due to the perpendicular component of the acceleration
GalDynPsr.pdotint.PdotShk(Ex_shk,P) # due to the Shklovskii effect
The total dynamically caused value of Pdot is the addition of the above three terms. The intrinsic value of the period derivative can be calculated by subtracting that sum from the measured value of the period derivative. GalDynPsr can perform this task for us by:
GalDynPsr.pdotint.Pdotint(Ex_pl,Ex_z,Ex_shk,Pdot_obs,P)
The total dynamical contribution in the period derivative, i.e. the sum of GalDynPsr.pdotint.PdotExpl(Ex_pl,P) and GalDynPsr.pdotint.PdotExz(Ex_z,P) can be printed as:
GalDynPsr.pdotint.PdotGal(Ex_pl,Ex_z,P)
One can calculate the uncertainity in the above period drivatives as:
GalDynPsr.pdotint.ErrPdotExpl(Ex_pl,errpl,P,sigP) # uncertainty in the planar contribution to the period derivative
GalDynPsr.pdotint.ErrPdotExz(Ex_z,errz,P,sigP) # uncertainty in the perpendicular contribution to the period derivative
GalDynPsr.pdotint.ErrPdotGal(Ex_pl,errpl,Ex_z,errz,P,sigP) # uncertainty in the sum of planar and perpendicular contributions to the period derivative
GalDynPsr.pdotint.ErrPdotShk(Ex_shk,errshk,P,sigP) # uncertainty in the Shklovskii term contribution to the period derivative
GalDynPsr.pdotint.ErrPdotint(Ex_pl,errpl,Ex_z,errz,Ex_shk,errshk,Pdot_obs,sigPdot_obs,P,sigP) # uncertainty in the intrinsic period derivative
* Remember the fact that GalDynPsr reads the value of the period in seconds. Orbital periods are usually measured in days (or sometimes in hours). One needs to convert those values in seconds before calling GalDynPsr modules.
### All the above points will be clearer from the following examples ################
##### Instructions to use GalDynPsr #####
# Call GalDynPsr as:-
import GalDynPsr
## Provide the following values in your code ###
# ldeg = Galactic logitude in degrees, bdeg = Galactic latitude in degrees, dkpc = distance to pulsar in kpc, sigl = error in ldeg, sigb = error in bdeg, sigd = error in dkpc
########### Example-1 #################
ldeg =20.0
sigl = 2.0
bdeg = 45.0
sigb = 0.2
dkpc = 1.2
sigd = 0.3
############# Extract important parameters say values of Rp (in kpc) and z (in kpc) and errors in them ##########
Rpkpc = GalDynPsr.read_parameters.Rpkpc(ldeg, sigl, bdeg, sigb, dkpc, sigd)
zkpc = GalDynPsr.read_parameters.z(ldeg, sigl, bdeg, sigb, dkpc, sigd)
ErrRp = GalDynPsr.read_parameters.ErrRp(ldeg, sigl, bdeg, sigb, dkpc, sigd)
Errz = GalDynPsr.read_parameters.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd)
################# Compute excess Shklovskii term. Here Exshk() calculates the dynamical contribution (or the excess term) due to Shklovskii effect and errShk() calculates the error in that #################################
#### We need to provide the values of proper motion in right ascesion direction (with the corresponding error) and proper motion in the declination direction (with the corresponding error)########
##mu_alpha = proper motion in Right Ascesion direction (mas/yr), mu_delta = proper motion in Declination direction (mas/yr), sigmua = error in mu_alpha, and sigmud = error in mu_delta
mu_alpha = 22.0
sigmua = 0.2
mu_delta = 1.2
sigmud = 0.3
ExcessSh = GalDynPsr.Shk.Exshk(dkpc, sigd, mu_alpha, sigmua, mu_delta, sigmud)
sigmaExcessSh = GalDynPsr.Shk.errShk(dkpc, sigd, mu_alpha, sigmua, mu_delta, sigmud)
### Example-1.1 ###
### Model A (can propagate error)######
Apl = GalDynPsr.modelA.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Az = GalDynPsr.modelA.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
Aplsigma = GalDynPsr.modelA.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the planar contribution to the excess term
Azsigma = GalDynPsr.modelA.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalA = Apl + Az
SigmatotalA = math.sqrt(Aplsigma*Aplsigma+Azsigma*Azsigma) # assuming no correlation between excepp_pl and excess_z
### Example-1.2 ###
### Model B (can propagate error)######
Bpl = GalDynPsr.modelB.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Bz = GalDynPsr.modelB.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
Bplsigma = GalDynPsr.modelB.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the planar contribution to the excess term
Bzsigma = GalDynPsr.modelB.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalB = Bpl + Bz
SigmatotalB = math.sqrt(Bplsigma*Bplsigma+Bzsigma*Bzsigma) # assuming no correlation between excepp_pl and excess_z
### Example-1.3 ###
### Model C (can propagate error)######
Cpl = GalDynPsr.modelC.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Cz = GalDynPsr.modelC.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
Cplsigma = GalDynPsr.modelC.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the planar contribution to the excess term
Czsigma = GalDynPsr.modelC.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalC = Cpl + Cz
SigmatotalC = math.sqrt(Cplsigma*Cplsigma+Czsigma*Czsigma) # assuming no correlation between excepp_pl and excess_z
### Example-1.4 ###
### Model D (can propagate error)######
Dpl = GalDynPsr.modelD.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Dz = GalDynPsr.modelD.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
Dplsigma = GalDynPsr.modelD.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the planar contribution to the excess term
Dzsigma = GalDynPsr.modelD.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalD = Dpl + Dz
SigmatotalD = math.sqrt(Dplsigma*Dplsigma+Dzsigma*Dzsigma) # assuming no correlation between excepp_pl and excess_z
### Example-1.5 ###
############ Model Ea (error only in parallel component) ###############
Eapl = GalDynPsr.modelEa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
SigmaEapl = GalDynPsr.modelEa.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the planar contribution to the excess term
Eaz = GalDynPsr.modelEa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
totalEa = Eapl+Eaz
### Example-1.6 ###
############ Model Eb (error only in parallel component) ###############
Ebpl = GalDynPsr.modelEb.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
SigmaEbpl = GalDynPsr.modelEb.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the planar contribution to the excess term
Ebz = GalDynPsr.modelEb.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
totalEb = Ebpl+Ebz
### Example-1.7 ###
############ Model Fa (error only in parallel component) ###############
Fapl = GalDynPsr.modelFa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
SigmaFapl = GalDynPsr.modelFa.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the planar contribution to the excess term
Faz = GalDynPsr.modelFa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
totalFa = Fapl+Faz
### Example-1.8 ###
############ Model Fb (error only in parallel component) ###############
Fbpl = GalDynPsr.modelFb.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
SigmaFbpl = GalDynPsr.modelFb.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the planar contribution to the excess term
Fbz = GalDynPsr.modelFb.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
totalFb = Fbpl+Fbz
### Example-1.9 ###
############ Model Ga (error in perpendicular component) ###############
Gapl = GalDynPsr.modelGa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Gaz = GalDynPsr.modelGa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
SigmaGaz = GalDynPsr.modelGa.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalGa = Gapl+Gaz
### Example-1.10 ###
############ Model Gb (error in perpendicular component) ###############
Gbpl = GalDynPsr.modelGb.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Gbz = GalDynPsr.modelGb.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
SigmaGbz = GalDynPsr.modelGb.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalGb = Gbpl+Gbz
### Example-1.11 ###
############ Model Ha (error in perpendicular component) ###############
Hapl = GalDynPsr.modelHa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Haz = GalDynPsr.modelHa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
SigmaHaz = GalDynPsr.modelHa.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalHa = Hapl+Haz
### Example-1.12 ###
############ Model Hb (error in perpendicular component) ###############
Hbpl = GalDynPsr.modelHb.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Hbz = GalDynPsr.modelHb.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
SigmaHbz = GalDynPsr.modelHb.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalHb = Hbpl+Hbz
### Example-1.13 ###
############ Model Ia (error in perpendicular component) ###############
Iapl = GalDynPsr.modelIa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Iaz = GalDynPsr.modelIa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
SigmaIaz = GalDynPsr.modelIa.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalIa = Iapl+Iaz
### Example-1.14 ###
############ Model Ib (error in perpendicular component) ###############
Ibpl = GalDynPsr.modelIb.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Ibz = GalDynPsr.modelIb.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
SigmaIbz = GalDynPsr.modelIb.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalIb = Ibpl+Ibz
### Example-1.15 ###
############ Model Ja (error in perpendicular component) ###############
Japl = GalDynPsr.modelJa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Jaz = GalDynPsr.modelJa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
SigmaJaz = GalDynPsr.modelJa.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalJa = Japl+Jaz
### Example-1.16 ###
############ Model Jb (error in perpendicular component) ###############
Jbpl = GalDynPsr.modelJb.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Jbz = GalDynPsr.modelJb.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
SigmaJbz = GalDynPsr.modelJb.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the error in the perpendicular contribution to the excess term
totalJb = Jbpl+Jbz
### Example-1.17 ###
############ Model Ka (full galpy, no error) ###############
Kapl = GalDynPsr.modelLa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Kaz = GalDynPsr.modelLa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
totalKa = Kapl+Kaz
### Example-1.18 ###
############ Model Kb (full galpy, no error) ###############
Kbpl = GalDynPsr.modelLa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Kbz = GalDynPsr.modelLa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
totalKb = Kbpl+Kbz
### Example-1.19 ###
############ Model La (full galpy, no error) ###############
Lapl = GalDynPsr.modelLa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Laz = GalDynPsr.modelLa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
totalLa = Lapl+Laz
### Example-1.20 ###
############ Model Lb (full galpy, no error) ###############
Lbpl = GalDynPsr.modelLa.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the planar contribution to the excess term
Lbz = GalDynPsr.modelLa.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd) #calculates the perpendicular contribution to the excess term
totalLb = Lbpl+Lbz
### Example-2 error estimation using Model La (full galpy, no error) ###
import numpy as np
import statistics as stat
from astropy import units as u
from astropy.coordinates import SkyCoord
import GalDynPsr
Number=50000 # number of instances to simulate
mul, sigmal = # define mean and standard deviation of l_degree
mub, sigmab = # define mean and standard deviation of b_degree
mud, sigmad = # define mean and standard deviation of d_kpc
# now simulate 'Number' instances of l, b and d
lsim = np.random.normal(mul, sigmal, Number)
bsim = np.random.normal(mub, sigmab, Number)
dsim = np.random.normal(mud, sigmad, Number)
# declare arrays for excess_parallel, excess_z and excess_total
Lapl = []
Laz = []
totalLa = []
# now run GalDynPsr for 'Number' times
count = 0
while (count < Number):
# randomly select l, b, d from the simulated values
x=np.random.random_integers(low=1, high=Number-1, size=3)
l=lsim[x[0]]
b=bsim[x[1]]
d=dsim[x[2]]
# for each value of l, b, d, calculate excess_parallel, excess_z and excess_total
Lapl.append(GalDynPsr.modelLa.Expl(l, sigmal, b, sigmab, d, sigmad))
Laz.append(GalDynPsr.modelLa.Exz(l, sigmal, b, sigmab, d, sigmad))
totalLa=np.add(Lapl,Laz)
count = count + 1
# now calculate mean and standard deviation of excess_parallel, excess_z and excess_total
print(stat.mean(Lapl), stat.stdev(Lapl), count)
print(stat.mean(Laz), stat.stdev(Laz), count)
print(stat.mean(totalLa), stat.stdev(totalLa), count)
#############################
######### For a Globular Cluster Pulsar #############################
#Check names of Globular Clusters from GC_Harris.txt and write the name of the cluster. We do not include errors in this case.
### Example 3 ###
#Getting Excess_pl for model A for the GC 47Tuc:-
#GalDynPsr.ReadGC.gcinput("GC Name") finds the GC in the GC_Harris list and calculates the corresponding galactic coordinates and distance to the GC. Input the arguments to any model as:- (GalDynPsr.ReadGC.ldeg,0.,GalDynPsr.ReadGC.bdeg,0.,GalDynPsr.ReadGC.dkpc,0.)
GalDynPsr.ReadGC.gcinput("47Tuc")
Apl = GalDynPsr.modelA.Expl(GalDynPsr.ReadGC.ldeg,0.,GalDynPsr.ReadGC.bdeg,0.,GalDynPsr.ReadGC.dkpc,0.)
Az = GalDynPsr.modelA.Exz(GalDynPsr.ReadGC.ldeg,0.,GalDynPsr.ReadGC.bdeg,0.,GalDynPsr.ReadGC.dkpc,0.)
# Similarly for other models.
#############################
######################### For Intrinsic period derivative calculations ###############################
## Parameters required for spin period derivative calcuations:- P = period in seconds, sigP = error in P in second, Pdot_obs = measured period derivative in ss^-1, sigPdot_obs = error in Pdot_obs
#PdotExpl(Excess_pl,P) calculates the planar contribution to the period derivative
#PdotExz(Excess_z,P) calculates the perpendicular contribution to the period derivative
#PdotGal(Excess_pl,Excess_z,P) calculates the sum of planar and perpendicular contributions to the period derivative
#PdotShk(Excess_shk,P) calculates the Shklovskii contribution to the period derivative
#Pdotint(Excess_pl,Excess_z,Excess_shk,Pdot_obs,P) calculates the intrinsic period derivative
#ErrPdotExpl(Excess_pl,error_Excess_pl,P,sigP) calculates the error in the planar contribution to the period derivative
#ErrPdotExz(Excess_z,error_Excess_z,P,sigP) calculates the error in the perpendicular contribution to the period derivative
#ErrPdotGal(Excess_pl,error_Excess_pl,Excess_z,error_Excess_z,P,sigP) calculates the error in the sum of planar and perpendicular contributions to the period derivative
#ErrPdotShk(Excess_shk,error_Excess_shk,P,sigP) calculates the error in the Shklovskii contribution to the spin period derivative
#ErrPdotint(Excess_pl,error_Excess_pl,Excess_z,error_Excess_z,Excess_shk,error_Excess_shk,Pdot_obs,sigPdot_obs,P,sigP) calculates the error in the intrinsic period derivative
### Example 4 ###
# Say using La
ldeg =20.0
sigl = 2.0
bdeg = 45.0
sigb = 0.2
dkpc = 1.2
sigd = 0.3
mu_alpha = 22.0
sigmua = 0.2
mu_delta = 1.2
sigmud = 0.3
# Example spin period and itsderivatives:-
#P = 0.003061844088094608
#sigP = 0.000000000000000015
#Pdot_obs = 0.959013e-20
#sigPdot_obs = 0.000014e-20
#Example orbital period and its derivatives (convert the period and corresponding error in seconds):-
P = 1.198512575184*24.0*3600.0
sigP = 0.000000000013*24.0*3600.0
Pdot_obs = 4.8e-14
sigPdot_obs = 1.1e-14
##
Ex_pl = GalDynPsr.modelA.Expl(ldeg, sigl, bdeg, sigb, dkpc, sigd)
Ex_z = GalDynPsr.modelA.Exz(ldeg, sigl, bdeg, sigb, dkpc, sigd)
errpl = GalDynPsr.modelA.Errpl(ldeg, sigl, bdeg, sigb, dkpc, sigd)
errz = GalDynPsr.modelA.Errz(ldeg, sigl, bdeg, sigb, dkpc, sigd)
Ex_shk = GalDynPsr.Shk.Exshk(dkpc, sigd, mu_alpha, sigmua, mu_delta, sigmud)
errshk = GalDynPsr.Shk.errShk(dkpc, sigd, mu_alpha, sigmua, mu_delta, sigmud)
PdotExpl = GalDynPsr.pdotint.PdotExpl(Ex_pl,P)
PdotExz = GalDynPsr.pdotint.PdotExz(Ex_z,P)
PdotGal = GalDynPsr.pdotint.PdotGal(Ex_pl,Ex_z,P)
PdotShk = GalDynPsr.pdotint.PdotShk(Ex_shk,P)
Pdotint = GalDynPsr.pdotint.Pdotint(Ex_pl,Ex_z,Ex_shk,Pdot_obs,P)
#Errors in spin period derivatives:-
ErrPdotExpl = GalDynPsr.pdotint.ErrPdotExpl(Ex_pl,errpl,P,sigP)
ErrPdotExz = GalDynPsr.pdotint.ErrPdotExz(Ex_z,errz,P,sigP)
ErrPdotGal = GalDynPsr.pdotint.ErrPdotGal(Ex_pl,errpl,Ex_z,errz,P,sigP)
ErrPdotShk = GalDynPsr.pdotint.ErrPdotShk(Ex_shk,errshk,P,sigP)
ErrPdotint = GalDynPsr.pdotint.ErrPdotint(Ex_pl,errpl,Ex_z,errz,Ex_shk,errshk,Pdot_obs,sigPdot_obs,P,sigP)
#==========================================================================================
##### Contents of the Package ####
Datafiles:
parameters.in: Input file contains values of different constants which are subject to change with improvements in measurements. User can change the values of constants if the need be.
GC_Harris.txt: Names and Galactic coordinates of the Globular Clusters.
README.txt: Contents of this README.md file inside package along with code files.
Description of different codes:
read_parameters.py: contains parameters (constants) used in the package
ExcessZ.py: Calculates az/c, both Nice-Taylor (1995) expression as well as our fit of HF2004
Excesspl.py: Calculates apl/c. Can take any value of the slope parameter dv/dR from the parameter file. b0 = (Vs/Rs)*(dV/dR) = 0 is the conventional.
modelA.py: Model A - Calculates apl_old/c (i.e. b0 = 0; DT91) and az_NT95.
modelB.py: Model B - Calculates apl_old/c (i.e. b0 = 0; DT91) and az_HF04.
modelC.py: Model C - Calculates apl_new/c (i.e. dv/dR = -0.2; Reid2014) and az_NT95.
modelD.py: Model D - Calculates apl_new/c (i.e. dv/dR = -0.2; Reid2014) and az_HF04.
modelEa.py: Model Ea - Calculates apl_old/c (i.e. b0 = 0; DT91) and ‘zforce’ in galpy (without BH)
modelEb.py: Model Eb - Calculates apl_old/c (i.e. b0 = 0; DT91) and ‘zforce’ in galpy (with BH)
modelFa.py: Model Fa - Calculates apl_new/c (i.e. dv/dR = -0.2; Reid2014) and ‘zforce’ in galpy (without BH)
modelFb.py: Model Fb - Calculates apl_new/c (i.e. dv/dR = -0.2; Reid2014) and ‘zforce’ in galpy (with BH)
modelGa.py: Model Ga - Calculates apl/c using Vp'/Vs from galpy (without BH) and az_NT95
modelGb.py: Model Gb - Calculates apl/c using Vp'/Vs from galpy (with BH) and az_NT95
modelHa.py: Model Ha - Calculates ‘Rforce’ in galpy (without BH) and az_NT95
modelHb.py: Model Hb - Calculates ‘Rforce’ in galpy (with BH) and az_NT95
modelIa.py: Model Ia - Calculates apl/c using Vp'/Vs from galpy (without BH) and az_HF04
modelIb.py: Model Ib - Calculates apl/c using Vp'/Vs from galpy (with BH) and az_HF04
modelJa.py: Model Ja - Calculates ‘Rforce’ in galpy (without BH) and az_HF04
modelJb.py: Model Jb - Calculates ‘Rforce’ in galpy (with BH) and az_HF04
modelKa.py: Model Ka - Calculates apl/c using Vp'/Vs from galpy (without BH) and ‘zforce’ in galpy (without BH)
modelKb.py: Model Kb - Calculates apl/c using Vp'/Vs from galpy (with BH) and ‘zforce’ in galpy (with BH)
modelLa.py: Model La - Calculates ‘Rforce’ in galpy (without BH) and ‘zforce’ in galpy (without BH)
modelLb.py: Model Lb - Calculates ‘Rforce’ in galpy (with BH) and ‘zforce’ in galpy (with BH)
Shk.py: Calculates the Shklovskii term d(mu_T*mu_T)/c. Takes user input of proper motion in right ascension and declination.
galpyMWpl.py: Using vcirc function in galpy to get Vp'(without BH) and eventually apl/c
galpyMWBHpl.py: Using vcirc function in galpy to get Vp'(with BH) and eventually apl/c
galpyMWRfo.py: Using evaluateRforces function in galpy(without BH) to get apl/c
galpyMWBHRfo.py: Using evaluateRforces function in galpy(with BH) to get apl/c
err_excesspl_Damour.py: Calculates error for Excess_pl (DT91)
err_excesspl_Reid.py: Calculates error for Excess_pl (Reid2014)
err_HFhigh.py: Calculates error for Excess_z (our fit for HF04 data)
err_HFlow.py: Calculates error for Excess_z (our fit for HF04 data)
err_NT.py: Calculates error for Excess_z (NT95)
err_Shklovskii.py: Calculates error for Excess_Shk
#==========================================================================================================================================
Though the description of models is provided in the paper, for the sake of completeness we are still providing that here.
Description of Models:-
Model-A: apl_old/c (i.e. b0 = 0; DT91) and Z_NT95
Model-B: apl_old/c (i.e. b0 = 0; DT91) and Z_HF04
Model-C: apl_new/c (i.e. dv/dR = -0.2; Reid2014) and Z_NT95
Model-D: apl_new/c (i.e. dv/dR = -0.2; Reid2014) and Z_HF04
Model-Ea: apl_old/c (i.e. b0 = 0; DT91) and ‘zforce’ in galpy (without BH)
Model-Eb: apl_old/c (i.e. b0 = 0; DT91) and ‘zforce’ in galpy (with BH)
Model-Fa: apl_new/c (i.e. dv/dR = -0.2; Reid2014) and ‘zforce’ in galpy (without BH)
Model-Fb: apl_new/c (i.e. dv/dR = -0.2; Reid2014) and ‘zforce’ in galpy (with BH)
Model-Ga: Using Vp'/Vs from galpy (without BH) and Z_NT95
Model-Gb: Using Vp'/Vs from galpy (with BH) and Z_NT95
Model-Ha: ‘Rforce’ in galpy (without BH) and Z_NT95
Model-Hb: ‘Rforce’ in galpy (with BH) and Z_NT95
Model-Ia: Using Vp'/Vs from galpy (without BH) and Z_HF04
Model-Ib: Using Vp'/Vs from galpy (with BH) and Z_HF04
Model-Ja: ‘Rforce’ in galpy (without BH) and Z_HF04
Model-Jb: ‘Rforce’ in galpy (with BH) and Z_HF04
Model-Ka: Using Vp'/Vs from galpy (without BH) and ‘zforce’ in galpy (without BH)
Model-Kb: Using Vp'/Vs from galpy (with BH) and ‘zforce’ in galpy (with BH)
Model-La: ‘Rforce’ in galpy (without BH) and ‘zforce’ in galpy (without BH)
Model-Lb: ‘Rforce’ in galpy (with BH) and ‘zforce’ in galpy (with BH)
Reid2014: Reid M. J., Menten K. M., Brunthaler A., Zheng X. W., Dame T. M., 2014, ApJ, 783, 130.
DT91: Damour T., Taylor J.H., 1991, ApJ, 366, 501.
NT95: Nice D.J., Taylor J.H., 1995, ApJ, 441, 429.
HF04: Holmberg J., Flynn C., 2004, MNRAS, 352, 440.
############################################################
نحوه نصب
نصب پکیج whl GalDynPsr-1.0.0:
pip install GalDynPsr-1.0.0.whl
نصب پکیج tar.gz GalDynPsr-1.0.0:
pip install GalDynPsr-1.0.0.tar.gz