How to select the psd-model when PEing GW150915 real data with PyCBC pycbc_inference executables?

There’re many optionable psd models for injecting a noise psd morphology in the PyCBC modules that is mainly produced from LALSimulation.
(like psd-model=aLGOxxxxxx as a parameter setting inside [data] section of ini. file)

From the LALSimulation document, I learned many noise psd or sensitivity models are the realizations of those from LIGO documents.
(and some of the documents seem unreachable for the unauthenticated)

So on earth which psd model should I choose to be the proper one for GW150914, then finally can I obtain the same parameter results as in the published GW150914 paper from LIGO official?

Or even should I trace back to the day of 150914 the status in which the LIGO detectors were working from the website provided by GWOSC?

Hi,

Great question! In general, for a specific event, it is best to use data around the event time to calculate a PSD.

If you are trying to replicate LVK results, I recommend using the PSD that is included in the PE Sample files. You can see an example of how to do that in this Jupyter Notebook

psd = data.psd["C01:IMRPhenomXPHM"]
fig = psd.plot(fmin=20)
ax = fig.gca()
ax.set_ylim(1e-48, 1e-40)

You should be able to access PSDs from the GW150914 PE samples that are in the GWTC-2.1 data release

If you like, you can also see how this is done in the Event Viewer App (though, I think the Jupyter notebook is easier to follow).

https://peviewer.igwn.org/?event1=GW150914

1 Like

Dear Jonah:

I have leafed through those links provided and they are really amazingly helpful to me !
Extremely Grateful to you, Mr Jonah !!

Recently, I’m trying to replicate the results of GW150914 from the data of LIGO, and actually got stuck.

Based on all my tryings, I got three questions:

  1. Should I estimate D_L, orbital plane inclination \theta_{JN} and sky loaction angular parameters \alpha and \delta at the same time together and with others? (For the degeneracy between the former two, and the indirect long way of gradually enhancing the precision of the sky-map credible region of the latter two)
  2. Should I marginalize anything for the case of GW150914? Maye this question is rooted from my still being ignorant about the purpose of marginalization.
  3. Should I estimate polarization \psi now for GW150914? As I noticed some sentences like this in a paper of event co-observed by LVK that “Finally we are able to estimate the polarization, for the joining of a third ground-borne detector…”

Would it be possible for you to give me some hints, Mr. Jonah?

Great thanks to what you have provided again, sir !
I will carefully check them for perhaps somewhere be the answers to my questions.


Background


I tried two examples with reference codes from PyCBC documents website:

  • One uses emcee_pt sampler, which means a standard MCMC sampler, and mainly shows the .ini files settings on the event GW150914.
  • Another one uses dynesty. which is for nested sampling, shows the settings for a time-marginalized model accompanied by also ra, dec, polarization, distance, phase(IMRPhenomD) marginalization.

Results:

  • emcee_pt sampler: Two components masses in detector frame seems reasonable, while other parameters are not satisfying. Especially the Luminosity Distance is estimated as around 600(±200).
  • dynesty sampler: the same as the picture at the page bottom of [PyCBC] Marginalized time model: Example with GW150914, which seems being the truth.

Main settings (not all)


dynesty case settings:

[model]
name = marginalized_time
low-frequency-cutoff = 30.0
marginalize_vector_params = tc, ra, dec, polarization
marginalize_phase = True
marginalize_distance = True



emcee_pt case settings:

ABOUT data.ini:
I used the straing data from LIGO’s H1 and L1 detectors:

[data]
frame-files = H1:H-H1_GWOSC_16KHZ_R1-1126257415-4096.gwf L1:L-L1_GWOSC_16KHZ_R1-1126257415-4096.gwf
channel-name = H1:GWOSC-16KHZ_R1_STRAIN L1:GWOSC-16KHZ_R1_STRAIN

This time segment was analysed:

analysis-start-time = -6
analysis-end-time = 2

Noise PSD was estimated 1024+6+8s ahead of -6s :

psd-estimation = median-mean
psd-start-time = -1038
psd-end-time = -6

I added 8s for data padding during the PSD estimation process (not sure if it should be treated like this):

psd-inverse-length = 8
psd-segment-length = 8
psd-segment-stride = 4
strain-high-pass = 15
pad-data = 8


ABOUT like.ini:
I chose the Gaussian-noise model

[model]
name = gaussian_noise
low-frequency-cutoff = 15.0

All parameters necessary for IMRPhenomPv2 was included:

[variable_params]
delta_tc =
coa_phase =
mass1 =
mass2 =
spin1_a =
spin1_azimuthal =
spin1_polar =
spin2_a =
spin2_azimuthal =
spin2_polar =
distance =
ra =
dec =
inclination =
polarization =


ABOUT sampler.ini:
I tried the emcee_pt sampler:

[sampler]
name = emcee_pt
nwalkers = 200
ntemps = 20

With burn_in setting as:

[sampler-burn_in]
burn-in-test = nacl & max_posterior