Reports until 18:00, Tuesday 02 April 2019
H1 CAL
sheila.dwyer@LIGO.ORG - posted 18:00, Tuesday 02 April 2019 - last comment - 20:57, Friday 05 April 2019(48178)
changes to pcal correction and delay in sensing measurement processing function, fixed error in DARM filter model

J Kissel, S Dwyer

We spent some time looking at sensing.py today, trying to understand better the systematic error in our sensing function model used for fitting the sensing function measurement.

Pcal Correction:

We also noticed that there is something called OMCDelay (and IOP delay which is insidetotalSensingDigital ) applied in Pcal correction, at first glace it doesn't seem right to apply an OMCDelay here. The default value of both of these delays is 0, and we think that omcDelay was intended to be a delay from sending pCal from the end station to  CalCS model in the corner.  Jeff filed a ticket : 12638  Since both of these are set to 0 delay, they don't impact the systematic error in the sensing or actuation functions. 

DARM excitation delay:

We committed this version of sensing.py at 7121.

DARM filters:

Non-image files attached to this report
Comments related to this report
ling.sun@LIGO.ORG - 10:14, Thursday 04 April 2019 (48240)

Lilli S, Jenne D, Jeff K,

I found the alog documenting when the DARM1 FM3 (resG) filter was turned on (47235). Jenne checked the exact time when the filter was turned on.

2019 Mar 01 23:38:37 UTC -- FM3 was turned on for the first time

2019 Mar 02 03:07-04:02 UTC -- It seems that during this period the filter was turned off again for some calibration testing.

I noticed that there was a measurement template during this time period: /ligo/svncommon/CalSVN/aligocalibration/trunk/Runs/O3/H1/Measurements/FullIFOSensingTFs/2019-03-01_H1_PCAL2DARM_TF_5t1100Hz_8min_AFTER_CHANGE.xml (2019-03-02 03:26:24 UTC)

It is not clear whether the model file modelparams_H1_20190301.py was generated using the measurements when FM3 was on or off. But for all the model files generated after Mar 01, I have added FM3 to the DARM1 filter bank.

sheila.dwyer@LIGO.ORG - 11:17, Friday 05 April 2019 (48265)

We have removed the excitation delay from the DARM EXC/ DARM IN2 measurement processing in sensing.py, because we think this was an error after Jenne looked into it more.  This will make our sensing function fit even worse, like the second attachment above.  

sensing.py is committed at revision 7143

jeffrey.kissel@LIGO.ORG - 20:57, Friday 05 April 2019 (48272)
J. Driggers, J. Kissel, L. Sun

TL;DR -- We also reverted the application of the full pcalCorrection mentioned above today, but now I think there's a bug in the way that the DC gain of the analog anti-aliasing filters (which is 0.99) are applied.

Here's why (as we were reminded by Evan and Shivaraj, but only just finally grokked today):

We take two measurements to obtain a complete measurement of the interferometer's response to differential arm displacement, C_meas,

  DARM_IN1       C_truth      (artifacts in PCAL_RXPD)
 ---------- = ------------- * -----------------------
 PCAL_RXPD    (1 + G)_truth   (artifacts in DARM loop)

  DARM_EXC                     
  --------  = (1 + G)_truth * (artifacts in DARM loop)
  DARM_IN2                  
which, as is indicated above, have some artifacts for which we need to correct.
When we multiply these together, the (artifacts in DARM loop) cancel, leaving

 DARM_IN1      DARM_EXC                            
 --------   *  --------  =  C_truth  *  (artifacts in PCAL_RXPD)
 PCAL_RXPD     DARM_IN2  

Now, we want to isolate C_truth. 

We know that (artifacts in PCAL_RXPD) are what we need to apply to turn the PCAL *channel* back in to DARM displacement,

(artifacts in PCAL_RXPD) = armSign_PCAL * (two poles at 1 Hz) * 1/(PCAL analog AA filter) * 1/(PCAL digital 65k-16k down sampling filter)

                         = armSign_PCAL * ("pcal dewhitening")* 1/(AAanalog_PCAL response * AAanalog_PCAL DC gain) * 1/(AAdigital_PCAL)

where AAanalog_PCAL has been measured to have a DC gain of 0.99 (import later), and AAdigital is a digital filter provided by the RCG with a DC gain of 1.00.

Thus, to get at C_meas alone, we take our two measured transfer functions and divide it by our known "pcalCorrection,"

 DARM_IN1      DARM_EXC                              armSign_PCAL * (two poles at 1 Hz)                       
 --------   *  --------  = C_truth * ----------------------------------------------------------------
 PCAL_RXPD     DARM_IN2             (AAanalog_PCAL response * AAanalog_PCAL DC gain) * AAdigital_PCAL 


           DARM_IN1      DARM_EXC      (AAanalog_PCAL response * AAanalog_PCAL DC gain) * AAdigial_PCAL                        
C_truth =  --------   *  --------  *   ----------------------------------------------------------------
           PCAL_RXPD     DARM_IN2                    armSign_PCAL *  (two poles at 1 Hz)   

We also have a model of this measurement of C_truth, which includes,

C_model =   (Optical Response = optical spring * coupled cavity pole) 
          * (Optical Gain) 
          * (IFO delay = (L/c) * single pole correction) 
          * (Uncompensated high frequency OMC DCPD Trans. Impedance Amp and Analog Whitening Filter Poles) 
          * (OMC DCPD analog AA filter response) * (OMC DCPD analog AA filter gain) 
          * ADC Gain 
          * (OMC DCPD digital 65k-16k down sampling filter)

where the (Optical Response) and (Optical Gain) are the only things we *don't* know or can't measure on the bench. We also *choose* to empirically measure the overall gain, from physical DARM displacement to DARM_IN1 counts, because it's easier to do that than measure the gains of all the electronics (i.e. the trans impedance amp, the whitening chassis, the analog AA, and other things between the "ideal" optical gain and the OMC DCPDs, i.e. between the power coming out of the SRC and the power on the DCPDs), because we can get very precise answer without uncertainty build-up from each term.
So, 
(overall gain) =   (Optical Gain [W/m]) 
                 * (loss between SRM to OMC DCPDs [W/W]) 
                 * (DCPD QE [A/W]) 
                 * (Trans. Impedance Gain [V/A]) 
                 * (Whitening Gain [V/V]) 
                 * (OMC analog AA Gain [V/V]) 
                 * (ADC gain [ct/ V]) 
                 * (whatever digital gains folks have implemented between the ADC channel and DARM_IN1 [ct/ct])

(and, of course, there are two DCPDs, so it's actually the split gain of the two paths recombined digitally in the real-time system, but you get the idea -- we want to fit for the gain of the whole path -- including the OMC DCPD's Analog AA filter gain.)

So we create a multi-parameter MCMC fit on C_truth to obtain the (overall gain) and (overall gain) but that means we must take out everything that we know, i.e.

C_MCMCinput = (Optical Response) * (overall gain) 
            = C_truth / (IFO delay * high freq. DCPD poles * AAanalog_DCPD response * AAdigital_DCPD)

Now here's where things get SNEAKY and CONFUSING.
     (1) We only have one model for the analogAA that we use for both the PCAL and the OMC DCPDs. Thus, analogAA_PCAL = analogAA_DCPD.
     (2) The RCG is doing identically the same thing for the digitalAA. Thus, digitalAA_PCAL = digitalAA_DCPD
That means that technically, we'd have to do *less* to create C_MCMCinput, because the (analogAA_PCAL * digitalAA_PCAL) from pcalCorrection cancel the (analogAA_DPCD * digitalAA_DCPD),

                                                         1
C_MCMCinput = C_truth * -----------------------------------------------------------------------------
                       (IFO delay * high freq. DCPD poles * AAanalog_DCPD response  * AAdigital_DCPD)  


               DARM_IN1      DARM_EXC      (AAanalog_PCAL response * AAanalog_PCAL DC gain) * AAdigial_PCAL                                      1
C_MCMCinput =  --------   *  --------  *   -----------------------------------------------------------------  *  -----------------------------------------------------------------------------
               PCAL_RXPD     DARM_IN2                       armSign_PCAL * (two poles at 1 Hz)                   (IFO delay * high freq. DCPD poles * AAanalog_DCPD response  * AAdigital_DCPD)    


               DARM_IN1      DARM_EXC                  (AAanalog_PCAL response * AAanalog_PCAL DC gain) * AAdigial_PCAL                                                  
C_MCMCinput =  --------   *  --------  *   -------------------------------------------------------------------------------------------------------------
               PCAL_RXPD     DARM_IN2      armSign_PCAL * (two poles at 1 Hz) * IFO delay * high freq. DCPD poles * AAanalog_DCPD response * AAdigital_DCPD  


               DARM_IN1      DARM_EXC                             AAanalog_PCAL DC gain              
C_MCMCinput =  --------   *  --------  *   ------------------------------------------------------------------------
               PCAL_RXPD     DARM_IN2       armSign_PCAL *  (two poles at 1 Hz) * IFO delay * high freq. DCPD poles

So -- you can get away with *only* dividing the measurement by *only* the frequency response (two poles at 1 Hz) and the (high freq. DCPD poles), which is what pyDARM was doing before, we just didn't realize it.. Note, I think pyDARM is still missing the IFO delay, but this works out to be 1 [us], so a small error on top of much bigger impactful things. By applying the full pcalCorrection, we were double counting the analogAA and digitalAA.

Note above, we're left with AAanalog_PCAL DC gain because we've chosen to fold the gain of the OMC DCPD analog AA into the fit of the overall gain.
I would have not have done this sneaky cancellation business, because there's confusion on the way -- BUT the math above gets you the true (overall gain) of the sensing function and the true (optical response).

BUT the DC Gain of the AAanalog_PCAL is not multiplied, it's divided, in pyDARM and I think it results in a missing multiplicative factor of AAanalog_PCAL DC gain = 0.99 in the construction of 1/C_GDS:
     (3) Again, the DC gain of *one* of analogAA models is *divided not multiplied* out of the what is input to the MCMC, and thus
     (4) The "optical gain" max aposteriori (MAP) that is spit out of the MCMC is then a fit to (overall gain) / (analogAA DC gain)^2.
AND
     (5) The user receives information to install 1/MAP into foton for the CAL-CS version of the 1/C correction, so this is wrong by (analogAA DC gain)^2 / (overall gain),
     (6) The sensing function model is reconstructed, it's done with the (overall gain) / (analogAA DC gain)^2 AND *one* of the analog models multiplied in *WITH ITS DC GAIN,* so we get model of C with (overall gain) / (analogAA DC gain), which means comparisons between (the MCMC input / model) are showing (overall gain) / (analogAA DC gain)^3,
     (7) gds corrections for the sensing function are *re*created ( from scratch without the model from (5) ), it's also created with *one* of the analog models *WITH ITS DC GAIN* multiplied, in so GDS corrections to C are 1 / (analogAA DC gain), and thus
     (8) 1/C_GDS = (analogAA DC gain) / C_truth

Let's look at 
    ^/trunk/Common/pyDARM/src/sensing.py   rev 7146

Inside the function "process_sensing_measurements,"

Lines 325-327      PCAL_tf = DARM_IN1 / PCAL_RXPD # essentially, later filtered for coherence, but overwritten with the same name
Lines 334-336       CLG_tf = DARM_EXC / DARM_IN1 # essentially, later filtered for coherence, but overwritten with the same name

Line 361           if modelAndDataPars[n].PCAL_arm == 'Y': sign = -1 # armSign_PCAL from above

                   # Here's where (3) gets applied
                   #                           armSign_PCAL       /        (two poles at 1 Hz)                           /               AAanalog_PCAL DC gain
Line 376           C_meas = PCAL_tf * CLG_tf * sign / signal.freqresp(sensProd.pcaldewhitening, 2.0*np.pi*PCAL_f)[1] / abs(signal.freqresp(sensProd.AAanalog, 2.0*np.pi*0.01)[1])
       
                   #              from above  /           high freq. DCPD poles
Line 380           opticalResponse = C_meas / signal.freqresp(sensProd.uncompensatedOMCDCPD, 2.0*np.pi*PCAL_f)[1]

and opticalResponse (assuming we don't correct for time dependence) is what's input into the MCMC fitting code because (in the same function), the named tuple "results" has its second slot filled with "opticalResponse",
Line 643           thisresult = results(PCAL_f, opticalResponse, unc, err, opticalResponse_mAperpm, err_mAperpm, opticalResponse_ref_tf)
Line 644           allresults.append(thisresult)

and then we zoom down to the function that runs the MCMC fitting, called sensing_MCMC, in which this "opticalResponse" is shoved into the MCMC hammer, 
Line 860           sampler = emcee.EnsembleSampler(nwalkers, ndim, lnprob_antispring, args=(PCAL_f, opticalResponse, err), threads=3)

from which the first of the maximum aposteriori values (MAP) values returned are the "optical gain," which we now know is (overall gain) / (analogAA DC gain)^2.
This MAP values is the "optical gain" that we stick in to any updated parameter file.

For CAL-CS, we just install the inverse of the "optical gain" MAP value, which we now know is (analogAA DC gain)^2 / (overall gain).

Then, this parameter file is used as input *back* in to sensing.py, where C_model is created that is eventually used to create filters for GDS. Here, the analogAA filter is multiplied back in *with its DC Gain*. See
Line 229          # FIXME the variable below is NAMED INCORRECTLY. It's actually C not 1/C. JSK 2019-04-04
Line 237          # gdscorrection = 1/ calcsResid2gdsInvSensOut = 1 / ( 1/C_pyDARM / 1/C_foton ) = (C_pyDARM / C_foton)
Line 238          calcsResid2gdsInvSensOut = signal.freqresp(serielZPK(serielZPK(uncompensatedOMCDCPD, AAanalog), omcDCPDtf), 2.0*np.pi*freq)[1] * AAanalog_delay   # analogAA is multiplied in WITH ITS DC GAIN
Line 239             * delay * singlePoleApproxDelayCorrection 
Line 240             * totalSensingDigital_freqresp 
Line 241             * 1.0/signal.dfreqresp(signal.ZerosPolesGain([],np.zeros(1),1,dt=1.0/2**14), 2.0*np.pi*freq/2**14)[1]  
Line 242             * pars.sensingSign 
Line 243             / IIRwarp_interp   

And this multiplied on to the data stream as 1/calcsResid2gdsInvSensOut (see res_corr_model on line 65 of GDS_FIR_filter_generation.py), so we have
     1                   1                  (analogAA DC gain)^2            1             (analogAA DC gain)
----------- * ------------------------  =   -------------------- *  ------------------ =  ------------------
 MAP value    calcsResid2gdsInvSensOut        (overall gain)        (analogAA DC gain)      (overall gain) 

Which means that at high frequency, GDS, is spitting out 0.99 / C_truth.

Also note that when the pcal correction is created, it's created *WITHOUT THE AA's DC GAIN,
    From above
    pcalCorrection = armSign_PCAL * ("pcal dewhitening")* 1/(AAanalog_PCAL response * AAanalog_PCAL DC gain) * 1/(AAdigital_PCAL)
and yet Line 208 (split for better labeling of parts) of sensing.py is
#                                     ("pcal dewhitening")           /               AAanalog_PCAL response                       /                 AAanalog_PCAL DC gainA           *         Adigital_PCAL        *  delay resolves to 0
pcalCorrection = signal.freqresp(pcaldewhitening, 2.0*np.pi*freq)[1] / (signal.freqresp(AAanalog, 2.0*np.pi*freq)[1]*AAanalog_delay/abs(signal.freqresp(AAanalog, 2.0*np.pi*0.01)[1]) * totalSensingDigital_freqresp * omcDelay_freqresp)


It's quite a kerfuffle. I've opened an FRS ticket for us to solve this issue -- see IIET Ticket 12642.