J. Kissel, L. Sun

We were all so proud to claim 1%/1deg systematic error in the PCAL 2 DELTAL EXTERNAL transfer function this morning -- see LHO aLOG 48040 -- but this required multiplying the actuator model in the front-end to be scaled by 0.95 (i.e. a 95% decrease). Turns out the problem lies in the detuned sensing function -- namely that a bug in the fitting code falsely reported no residual frequency-dependent systematic error between model and measurement, when -- in fact -- the model was bogus, and there was indeed a substantial residual -- to the tune of ~5% below 20 Hz -- consistent with the features we were trying to remove by scaling the actuator.

Since the assumption is that we've falsely put this 5% on the actuator, we are over reporting our sensitivity by 5% below 20 Hz. Since this is so low in frequency, it shouldn't have much impact on the BNS range. However, there must be good agreement between measurement and model (i.e. frequency dependent errors of all functions must be below ~1%) in order for the time dependent correction factor calculations to work, because the calibration lines at all frequencies (but especially those at 20 Hz) are multiplied by loop correction factors produced by the model.

I'll explain in more detail below, but the conclusion is that 
    - In order to obtain better agreement, we must better measure the interferometer detuning [namely, measure below 5 Hz], in order to better constrain the MCMC fitting of the spring frequency and Q.
    - Once we're happy with that, we must update the CAL-CS inverse sensing function filter (again).
    - Then we can re-verify that 
        - we don't need any fudged factor of 0.95 on the CAL-CS implementation of the actuator model,
        - the PCAL2DELTAL transfer function is flat, and
        - the CAL-CS implementation of the model matches the model itself with good fidelity

We have received the OK from Keita, the run manager to execute this plan tomorrow (Saturday, 2019-03-30) starting around 11:00a [assuming the IFO is up and running].


DETAILS

A few days ago, the primary sensing function code in the pyDARM infrastructure,
    /ligo/svncommon/CalSVN/aligocalibration/trunk/Common/pyDARM/src/sensing.py
was modified to be able to handle a pro-spring detuning which now both L1 and H1 see (e.g. LHO aLOG 47941 or LHO aLOG 42791; of course, L1 sees it to a lesser extent because H1's point absorbers). H1's spring is rather prominent in our current configuration given our choice of SR3 heating (see e.g. LHO aLOG 47604), only decided upon last week for improved noise.

This core code is one file with 1000+ lines composed of many functions including,
    (1) process the collection of sensing function measurements, 
    (2) uses an MCMC algorithm to fit that measurement, 
    (3) generates a model of the sensing function based on the most likely, maximum a posteriori (MAP), values of the MCMC fit, plus a foton design string to invert it
    (4) divides that new fit model by measurement to to produce an estimate of residual frequency-dependent systematic error, and
    (5) uses a Gaussian process algorithm to fit the remaining frequency-dependent systematic error to a posterior distribution of frequency-dependent functions.
It also has functions that are mere tools like,
    (A) create the optical plant frequency response of the now-standard detuned-FPDRMI given input of the cavity pole frequency, and detuned optical spring frequency and Q
    (B) assembling all the details of the full sensing function on to an input optical plant frequency response
    (C) creating an estimate of the uncertainty given an input coherence
    (D) creating a Gaussian probability density function from input data vector, with mean, and 1-sigma uncertainty
    (E) creating a optical plant seed function/equation over which to run the MCMC fit from input parameters
but it was (A) and (E) that were modified, and for-better-or-worse, the equations to create the optical response for these two functions were not common to both.

In order to invoke the new functionality model of a pro-spring (instead of an anti-spring which had only ever been seen at H1 during O1 and O2), one set the parameter for the spring frequency to be imaginary inside the collection of pyDARM model parameters, e.g.
    /ligo/svncommon/CalSVN/aligocalibration/trunk/Runs/O3/H1/params/modelparams_H1_20190328.py

However, in the modification of (A), the code checked if the input argument was imaginary, and then if so applied the following math for the input arguments detuneSpringFreq (now imaginary) and detuneSpringQ: 
    Line 62: detuneFunction = signal.TransferFunction([1,0,0],[1,2.0*np.pi*detuneSpringFreq/detuneSpringQ,-(2.0*np.pi*detuneSpringFreq)**2])

We discovered this afternoon, that with an imaginary input argument for the spring frequency, this math produces a bogus optical spring response only after comparing the CAL-CS implementation of a spring with fs = 5.620 Hz and Q = 5.201 -- which had been entered into foton directly from the MAP of the MCMC by humans (LHO aLOG 47941), bypassing the (currently broken for pro-spring) foton string from function (3) -- against the pyDARM produced optical plant.

This comparison is the first .pdf attachment, 2019-03-29_H1_C_pyDARM_vs_CALCS_beforemodifiedprospring.pdf.

For a comparison of what the correct math should be for a pro spring, see the second .pdf attachment, 2019-03-29_bogus_spring_demo.pdf.

This bug was fixed by changing the same line of code to be 
    Line 62: detuneFunction = signal.TransferFunction([1,0,0],[1,2.0*np.pi*abs(detuneSpringFreq)/detuneSpringQ,(2.0*np.pi*abs(detuneSpringFreq))**2])

With this fix, the pyDARM code produced a virtually identical response to foton, as expected
The new comparison is shown in the third .pdf attachment, 2019-03-29_H1_C_pyDARM_vs_CALCS_modifiedprospring.pdf

But -- what does this all mean for the systematic error in the sensing function?

Well, after quadruple checking, it turns out function that creates the MCMC fit seeding function (E), was split into two functions (one for pro- and one for anti- spring) and insensitive to this bug. The pro-spring function provided the correct pro-spring formula,
     Line 739: est_opt_response_tf = mu[0] / (1.0 + 1j*xdata/mu[1]) * (xdata**2/(xdata**2 - mu[2]**2 + 1*xdata*mu[2]*mu[3])) * np.exp( -2.0*np.pi*1j * mu[4]*1.0e-6 * xdata )
with an inputs of a frequency vector xdata and five-element vector of mean value priors, mu = [optical gain, cavity pole frequency, spring frequency, spring Q, and a delay], given that the values for the mu are hard-coded on line 786 (yuck!). The hard-coded prior for the spring frequency is a real, positive 5.0 Hz, i.e. expecting an anti-spring.

That meant that the updated MAP values for these parameters from the MCMC fit function, (2), were "right" but the function used to create a model response from those values, i.e. (A), was producing the wrong frequency response that does not accurately represent the MAP values.

*That* meant that when this MAP updated model was divided from the measurement, in (4), the amount of residual frequency-dependent systematic error was wrong, and coincidentally small.

To see this, look at the fourth and fifth .pdf attachments side-be-side which compare MCMC fit results against model,
    - before the bug fix 2019-03-29_H1_sensingFunction_beforemodifiedprospring.pdf
    - after the bug fix 2019-03-29_H1_sensingFunction_modifiedprospring.pdf

Since we now are using the fixed, committed function (A), the latter is now physically correct, but sadly, it (finally truthfully) reveals that the MCMC fit to the detuning is actually "quite" poor. In other words, there is not enough low-frequency data to properly constrain the MCMC fit.

The right panels of the 3rd page, or all of the fourth page show the true systematic error between model and measurement -- and this true systematic error is consistent with a 5% error [1.05] around 20 Hz, which is consistent with having to scale the actuator down by 5%.

You'll notice that in today's broadband and sweep PCAL2DARM transfer functions, well-below the sensor / actuator cross-over frequency, and especially below 20 Hz, the transfer function is asymptoting to 0.95 -- exactly the "correction" that I put in.

Using some toy math, you can see this immediately assuming that we don't model the sensing measurement, C_meas, well, then we could fake the answer at the cross-over frequency and above by adjusting the gain on A:
    dL_real = (1/C_meas) d_err + A_meas * d_ctrl

    dL_real = (1/C_model) (C_model/C_meas) * d_err + A_meas * d_ctrl
   
    (C_meas/C_model) dL_real = (1/C_model) * d_err + (C_meas/C_model) * A_meas * d_ctrl

*phew* What an insidious bug!!

So, again, we must measure *more* of the detuned spring shape (i.e. to lower frequency) in order to nail down the frequency and Q -- which will help nail down the optical gain, and get everything right in the end.

We'll execute the above plan tomorrow.