18 Wrapper class for the PulsarParameters structure, so that it can be accessed
19 in a dictionary-like way.
23 from __future__
import division, absolute_import, print_function
28 from six
import string_types
32 from astropy
import units
as u
43 "DM": u.pc / (u.cm) ** 3,
44 "DM1": u.pc / (u.cm**3 * u.yr),
66 "GLF2": u.Hz / u.s**2,
74 "EPS1": u.dimensionless_unscaled,
75 "EPS2": u.dimensionless_unscaled,
84 "SINI": u.dimensionless_unscaled,
87 "DR": u.dimensionless_unscaled,
88 "DTHETA": u.dimensionless_unscaled,
89 "SHAPMAX": u.dimensionless_unscaled,
101 "D_AOP": 1.0 / u.rad,
111 "H0": u.dimensionless_unscaled,
112 "APLUS": u.dimensionless_unscaled,
113 "ACROSS": u.dimensionless_unscaled,
116 "COSIOTA": u.dimensionless_unscaled,
118 "C22": u.dimensionless_unscaled,
119 "C21": u.dimensionless_unscaled,
122 "CGW": u.dimensionless_unscaled,
124 "COSTHETA": u.dimensionless_unscaled,
126 "I21": u.dimensionless_unscaled,
127 "I31": u.dimensionless_unscaled,
128 "Q22": u.kg * u.m**2,
129 "H0_F": u.dimensionless_unscaled,
130 "HPLUS": u.dimensionless_unscaled,
131 "HCROSS": u.dimensionless_unscaled,
134 "HSCALARB": u.dimensionless_unscaled,
135 "HSCALARL": u.dimensionless_unscaled,
138 "HVECTORX": u.dimensionless_unscaled,
139 "HVECTORY": u.dimensionless_unscaled,
142 "HPLUS_F": u.dimensionless_unscaled,
143 "HCROSS_F": u.dimensionless_unscaled,
144 "PSITENSOR_F": u.rad,
145 "PHI0TENSOR_F": u.rad,
146 "HSCALARB_F": u.dimensionless_unscaled,
147 "HSCALARL_F": u.dimensionless_unscaled,
148 "PSISCALAR_F": u.rad,
149 "PHI0SCALAR_F": u.rad,
150 "HVECTORX_F": u.dimensionless_unscaled,
151 "HVECTORY_F": u.dimensionless_unscaled,
152 "PSIVECTOR_F": u.rad,
153 "PHI0VECTOR_F": u.rad,
154 "TRANSIENTSTARTTIME": u.s,
166 "PMRA": u.mas / u.yr,
167 "PMDEC": u.mas / u.yr,
170 "PMELONG": u.mas / u.yr,
171 "PMELAT": u.mas / u.yr,
174 "PMBETA": u.mas / u.yr,
175 "PMLAMBDA": u.mas / u.yr,
185 "OMDOT": u.deg / u.yr,
194 "D_AOP": 1.0 / u.arcsec,
201 "TRANSIENTSTARTTIME": u.d,
216 A class to wrap the SWIG-wrapped lalpulsar.PulsarParameters structure.
218 This class lets you access the structure in a more Pythonic way, as well as providing
219 a nice format for holding pulsar (`.par`) parameter files.
221 The class can be used to set numerical values (double precision, unsigned integers), strings,
222 or vectors of floating point values, e.g.:
224 >>> from lalpulsar.PulsarParametersWrapper import PulsarParametersPy
225 >>> pppy = PulsarParametersPy() # an empty structure
226 >>> pppy['DECJ'] = 0.23 # set a numerical value
227 >>> pppy['BINARY'] = 'BT' # set a string value
228 >>> pppy['F'] = [10.2, 1.4e-11] # set a vector of float values
231 pp (PulsarParameters, str): a lalpulsar.PulsarParameters structure, or a string giving the
232 path to a TEMPO-style (`.par`) pulsar parameter file. If nothing is given then an empty
233 lalpulsar.PulsarParameters structure is created. The `read()` method can subsequently
234 be used to read in a `.par` file, or parameters can be added.
237 An example of initialising the class with a previously created `lalpulsar.PulsarParameters`
241 >>> from lalpulsar.PulsarParametersWrapper import PulsarParametersPy
242 >>> # read in a pulsar parameter file
243 >>> pp = lalpulsar.ReadTEMPOParFile('apulsar.par')
244 >>> # view as a PulsarParametersPy object
245 >>> pppy = PulsarParametersPy(pp)
247 The same thing could be achieved more directly using:
249 >>> pppy = PulsarParametersPy('apulsar.par')
251 or, equivalently with:
253 >>> pppy = PulsarParametersPy()
254 >>> pppy.read('apulsar.par')
260 _pulsarparameters =
None
268 if not isinstance(pp, lalpulsar.PulsarParameters)
and isinstance(
271 if os.path.isfile(pp):
275 raise ValueError(
"Input string does not point to a file")
276 elif isinstance(pp, lalpulsar.PulsarParameters):
280 "Expected 'lalpulsar.PulsarParameters' type, string, or None"
293 Get value from pulsar parameters
302 if key[-4:].upper() ==
"_ERR":
311 sname = re.sub(
r"_\d",
"", tkey)
if "_" in tkey
else re.sub(
r"\d",
"", tkey)
314 if sname != tkey
and tkey[0:2] !=
"GL":
318 int(tkey.split(
"_")[-1])
if "_" in tkey
else int(tkey[len(sname) :])
330 if not lalpulsar.PulsarCheckParam(self.
_pulsarparameters_pulsarparameters, tkey):
334 if not lalpulsar.PulsarCheckParam(self.
_pulsarparameters_pulsarparameters, tkey):
338 ptype = lalpulsar.PulsarGetParamType(self.
_pulsarparameters_pulsarparameters, tkey)
340 if ptype == lalpulsar.PULSARTYPE_REAL8_t:
342 value = lalpulsar.PulsarGetREAL8Param(self.
_pulsarparameters_pulsarparameters, tkey)
344 value = lalpulsar.PulsarGetREAL8ParamErr(self.
_pulsarparameters_pulsarparameters, tkey)
345 elif ptype == lalpulsar.PULSARTYPE_REAL8Vector_t:
348 tmpvalue = lalpulsar.PulsarGetREAL8VectorParam(
355 value = lalpulsar.PulsarGetREAL8VectorParamIndividual(
360 tmpvalue = lalpulsar.PulsarGetREAL8VectorParamErr(
363 value = tmpvalue.data
365 value = lalpulsar.PulsarGetREAL8VectorParamErrIndividual(
368 elif ptype == lalpulsar.PULSARTYPE_string_t:
370 value = lalpulsar.PulsarGetStringParam(self.
_pulsarparameters_pulsarparameters, tkey)
372 raise ValueError(
"String-type cannot have an error")
373 elif ptype == lalpulsar.PULSARTYPE_UINT4_t:
375 value = lalpulsar.PulsarGetUINT4Param(self.
_pulsarparameters_pulsarparameters, tkey)
377 raise ValueError(
"UINT4-type cannot have an error")
379 raise ValueError(
"Unrecognised type")
385 Set the value of a key
392 if isinstance(value, float):
393 lalpulsar.PulsarAddREAL8Param(self.
_pulsarparameters_pulsarparameters, key, value)
394 elif isinstance(value, string_types):
395 lalpulsar.PulsarAddStringParam(self.
_pulsarparameters_pulsarparameters, key, value)
396 elif isinstance(value, int):
398 lalpulsar.PulsarAddREAL8Param(self.
_pulsarparameters_pulsarparameters, key, float(value))
400 lalpulsar.PulsarAddUINT4Param(self.
_pulsarparameters_pulsarparameters, key, value)
401 elif isinstance(value, list)
or isinstance(value, np.ndarray):
402 tarray = lal.CreateREAL8Vector(len(value))
403 for i, tv
in enumerate(value):
404 if isinstance(tv, float):
407 raise ValueError(
"Non-float value in list or array")
408 lalpulsar.PulsarAddREAL8VectorParam(self.
_pulsarparameters_pulsarparameters, key, tarray)
410 raise ValueError(
"Data-type not one of know types")
414 Convert parameter values to equivalent dimensional versions.
417 name (str): the name of the parameter to convert
418 value: the value to the parameter to convert
421 :class:`astropy.unit.Unit`: a unit class with dimensions for a float parameters, or a
422 list containing unit classes for a list or :class:`numpy.ndarray`
425 if isinstance(value, string_types):
431 ppunit = u.dimensionless_unscaled
433 ppunit = PPUNITS[uname]
435 if isinstance(value, np.ndarray)
or isinstance(value, list):
439 cvalue.append(v * ppunit)
440 if uname
in [
"F",
"FB",
"P"]:
443 cvalue = value * ppunit
449 Convert from PulsarParameter units to TEMPO-par-file units
452 name (str): a parameter name
453 value (float, list, :class:`numpy.ndarray`): the value of the parameter
454 iserr (bool): state whether where converting the parameter's error value
457 :class:`astropy.unit.Unit`: a unit class with dimensions for a float parameters, or a
458 list containing unit classes for a list or :class:`numpy.ndarray`
461 from astropy.coordinates
import ICRS, Angle
462 from astropy.time
import Time
467 binaryunits = [
"XDOT",
"PBDOT",
"EPS1DOT",
"EPS2DOT",
"XPBDOT"]
483 "TRANSIENTSTARTTIME",
490 ppunit = PPUNITS[uname]
493 if uname
in TEMPOUNITS:
495 tempounit = TEMPOUNITS[uname]
497 if uname
not in TEMPOERRUNITS:
498 tempounit = TEMPOUNITS[uname]
500 tempounit = TEMPOERRUNITS[uname]
503 if uname
in binaryunits:
505 if abs(value) / 1e-12 > 1e-7:
508 if isinstance(value, string_types):
511 return value * u.dimensionless_unscaled
516 if ppunit == tempounit
or tempounit
is None:
519 if uname
in [
"RA",
"RAJ"]:
522 c = ICRS(pvalue, 0.0 * u.rad)
523 cvalue = c.ra.to_string(unit=tempounit, sep=
":", precision=12, pad=
True)
525 angle = Angle(pvalue)
527 angle.hms[0] * (60.0**2) + angle.hms[1] * 60.0 + angle.hms[2]
529 elif uname
in [
"DEC",
"DECJ"]
and not iserr:
530 c = ICRS(0.0 * u.rad, pvalue)
531 cvalue = c.dec.to_string(unit=tempounit, sep=
":", precision=12, pad=
True)
532 elif uname
in epochpars
and not iserr:
533 if isinstance(pvalue, list):
535 for i, pv
in enumerate(pvalue):
536 t = Time(pv, format=
"gps", scale=
"tt")
537 cvalue.append(t.mjd * tempounit)
539 t = Time(pvalue, format=
"gps", scale=
"tt")
540 cvalue = t.mjd * tempounit
541 elif uname
in binaryunits
and not iserr:
543 if abs(pvalue.value) / 1e-12 > 1e-7:
544 pvalue.value /= 1e-12
545 cvalue = pvalue.to(tempounit)
548 if isinstance(pvalue, list):
551 for i, pv
in enumerate(pvalue):
552 cvalue.append(pv.to(tempounit))
553 if uname
in [
"F",
"FB",
"P"]:
558 cvalue = pvalue.to(tempounit)
562 def parameter(self, name, withunits=False, tempounits=False):
564 Return the parameter given by name.
567 name (str): the name of the parameter to return
568 withunits (bool): if True return the parameter in a form with its appropriate units
569 tempounits (bool): of True return the parameter converted into the units required in a
570 TEMPO-style parameter file
578 if name.upper()[-4:] ==
"_ERR":
579 uname = name.upper()[:-4]
596 Return a list of the parameter names stored in the PulsarParameters structure
603 tname = thisitem.name
607 thisitem = thisitem.next
616 Return the values of each parameter in the structure
621 keys = self.
keyskeys()
625 ptype = lalpulsar.PulsarGetParamType(self.
_pulsarparameters_pulsarparameters, key)
627 if ptype == lalpulsar.PULSARTYPE_REAL8_t:
628 value = lalpulsar.PulsarGetREAL8Param(self.
_pulsarparameters_pulsarparameters, key)
629 elif ptype == lalpulsar.PULSARTYPE_REAL8Vector_t:
630 tmpvalue = lalpulsar.PulsarGetREAL8VectorParam(
636 elif ptype == lalpulsar.PULSARTYPE_string_t:
637 value = lalpulsar.PulsarGetStringParam(self.
_pulsarparameters_pulsarparameters, key)
638 elif ptype == lalpulsar.PULSARTYPE_UINT4_t:
639 value = lalpulsar.PulsarGetUINT4Param(self.
_pulsarparameters_pulsarparameters, key)
641 raise ValueError(
"UINT4-type cannot have an error")
643 raise ValueError(
"Could not find {} in strcuture".format(key))
645 tvalues.append(value)
651 Return the contents (not error at the moment) of the structure as a dictionary
656 for tpair
in self.
itemsitems():
657 tdict[tpair[0]] = tpair[1]
663 Return list of item tuples for each parameter in the structure
666 tkeys = self.
keyskeys()
667 tvalues = self.
valuesvalues()
671 for tk, tv
in zip(tkeys, tvalues):
672 titems.append((tk, tv))
676 def read(self, filename):
678 Read a TEMPO-style parameter file into a PulsarParameters structure
681 filename (str): the path to the pulsar `.par` file.
688 pp = lalpulsar.ReadTEMPOParFile(filename)
692 "Problem reading in pulsar parameter file '{}'".format(filename)
699 Return the error value for a particular parameter
702 name (str): the name of the parameter
706 if name.upper()[-4:] ==
"_ERR":
707 return self[name.upper()]
709 uname = name.upper() +
"_ERR"
716 Return the "fit flag" (a 1 or 0 depending whether the parameter with fit for by TEMPO(2)
720 name (str): the name of the parameter
723 if name.upper()
not in self.
keyskeys():
726 fitflag = lalpulsar.PulsarGetParamFitFlagAsVector(self.
_pulsarparameters_pulsarparameters, name)
728 if len(fitflag.data) > 1
or isinstance(self[name.upper()], (list, np.ndarray)):
731 return fitflag.data[0]
735 Convert the PulsarParameter structure to a string in the format of a
736 TEMPO-style `.par` file.
739 precision (int): the number of decimal places for an output value
742 str: the contents in TEMPO-format
747 max([len(kn)
for kn
in self.
keyskeys()]) + 5
750 outputstr =
"{{name: <{0}}}{{value: <{1}}}{{fitflag}}\t{{error}}".format(
759 for item
in self.
itemsitems():
763 if key
in [
"WAVESIN",
"WAVECOS"]
and parsedwaves:
771 if key
in TEMPOUNITS:
773 if evalue
is not None:
777 if isinstance(uvalue, list):
780 for tv, te
in zip(uvalue, uevalue):
781 tvalue.append(tv.value)
782 tevalue.append(te.value)
783 elif isinstance(uvalue, string_types):
785 tevalue = uevalue.value
787 tvalue = uvalue.value
788 tevalue = uevalue.value
794 if evalue
is not None:
799 if isinstance(tvalue, list)
or isinstance(tvalue, np.ndarray):
820 if key
in [
"WAVESIN",
"WAVECOS"]:
826 wavecos = self[
"WAVECOS"]
828 wavesin = self[
"WAVESIN"]
831 for ws, wc
in zip(wavesin, wavecos):
832 precstrs =
"{{0:.{}f}}".format(precision)
833 if ws < 1e-6
or ws > 1e6:
835 precstrs =
"{{0:.{}e}}".format(precision)
837 precstrc =
"{{0:.{}f}}".format(precision)
838 if wc < 1e-6
or wc > 1e6:
840 precstrc =
"{{0:.{}e}}".format(precision)
843 outputdic[
"name"] =
"WAVE{}{}".format(idxsep, idxoffset)
845 outputdic[
"value"] =
"{}\t{}".format(
846 precstrs.format(ws), precstrc.format(wc)
848 outputdic[
"fitflag"] =
""
849 outputdic[
"error"] =
""
851 parstr += outputstr.format(**outputdic).strip() +
"\n"
853 for tv, te, tf
in zip(tvalue, tevalue, fitflag):
854 precstr =
"{{0:.{}f}}".format(precision)
855 if tv < 1e-6
or tv > 1e6:
857 precstr =
"{{0:.{}e}}".format(precision)
859 precstre =
"{{0:.{}f}}".format(precision)
860 if te < 1e-6
or te > 1e6:
862 precstre =
"{{0:.{}e}}".format(precision)
865 outputdic[
"name"] =
"{}{}{}".format(key, idxsep, idxoffset)
867 outputdic[
"value"] = precstr.format(tv)
868 outputdic[
"fitflag"] =
"1" if tf == 1
else ""
869 outputdic[
"error"] = precstre.format(te)
if te != 0.0
else ""
871 parstr += outputstr.format(**outputdic).strip() +
"\n"
873 if isinstance(tvalue, float)
or key
in [
"RA",
"RAJ",
"DEC",
"DECJ"]:
874 if isinstance(tvalue, float):
875 precstr =
"{{0:.{}f}}".format(precision)
876 if tvalue < 1e-6
or tvalue > 1e6:
878 precstr =
"{{0:.{}e}}".format(precision)
880 ovalue = precstr.format(tvalue)
884 precstre =
"{{0:.{}f}}".format(precision)
886 if tevalue
is not None:
888 if tevalue < 1e-6
or tevalue > 1e6:
890 precstre =
"{{0:.{}e}}".format(precision)
892 oevalue = precstre.format(tevalue)
894 if fitflag
is not None:
901 outputdic[
"name"] = key
902 outputdic[
"value"] = ovalue
903 outputdic[
"fitflag"] = ofitflag
904 outputdic[
"error"] = oevalue
906 parstr += outputstr.format(**outputdic).strip() +
"\n"
910 def pp_to_par(self, filename, precision=19):
912 Output the PulsarParameter structure to a `.par` file.
915 filename (str): the path to the output file
916 precision (int): the number of decimal places for an output value
920 fp = open(filename,
"w")
922 raise IOError(
"Could not open file '{}' for writing".format(filename))
924 fp.write(self.
pp_to_strpp_to_str(precision=precision))
929 Return the PulsarParameters structure
936 Create a copy of the parameters.
940 memo[
id(self)] = newpar
941 lalpulsar.PulsarCopyParams(self.
PulsarParametersPulsarParameters(), newpar.PulsarParameters())
A class to wrap the SWIG-wrapped lalpulsar.PulsarParameters structure.
def __init__(self, pp=None)
def read(self, filename)
Read a TEMPO-style parameter file into a PulsarParameters structure.
def pp_to_par(self, filename, precision=19)
Output the PulsarParameter structure to a .par file.
def get_fitflag(self, name)
Return the "fit flag" (a 1 or 0 depending whether the parameter with fit for by TEMPO(2) in the ....
def convert_to_tempo_units(self, name, value, iserr=False)
Convert from PulsarParameter units to TEMPO-par-file units.
def __deepcopy__(self, memo)
Create a copy of the parameters.
def __setitem__(self, key, value)
Set the value of a key.
def __getitem__(self, key)
Get value from pulsar parameters.
def parameter(self, name, withunits=False, tempounits=False)
Return the parameter given by name.
def pp_to_str(self, precision=19)
Convert the PulsarParameter structure to a string in the format of a TEMPO-style ....
def as_dict(self)
Return the contents (not error at the moment) of the structure as a dictionary.
def keys(self)
Return a list of the parameter names stored in the PulsarParameters structure.
def PulsarParameters(self)
Return the PulsarParameters structure.
def get_error(self, name)
Return the error value for a particular parameter.
def convert_to_units(self, name, value)
Convert parameter values to equivalent dimensional versions.
def items(self)
Return list of item tuples for each parameter in the structure.
def values(self)
Return the values of each parameter in the structure.