I did a raw calibration of MCL and GUR. Accuracy is a factor of 2.

GUR path : 800 V/m/s => readout box (G~100) => ADC (0.7 mV/count)

MCL path : laser 1 MHz / V, cavity length ~ 25 m

I measured feedback signal before the laser with SR and avoided whitening filters for MC_F.

 PMC transmission is oscillating  in the range 0.5 - 0.85. PMC PZT voltage is 1-2 V.

FSS slow controls was -2.5 V. I adjusted it to 0 and PMC stabilized. PMC PZT voltage is 128, transmission is 0.845.

But most probably, slow control will drift again.

As I've spent many hours trying to determine the error in my C code for online filter I decided to write about it to prevent people from doing it again.

I have a C function that was tested offline. I compiled and installed it on the front end machine without any errors. When I've restarted the model, it did not run.

I modified the function the following way

I've adjusted input parameters such that STATEMENT was always true. However the model either started or not depending on the code after if statement. It turned out that the model could not start because of the following lines

cosine[1] = 1.0 - 0.5*a*a + a*a*a*a/24 - a*a*a*a*a*a/720 + a*a*a*a*a*a*a*a/40320;
sine[1] = a - a*a*a/6 + a*a*a*a*a/120 - a*a*a*a*a*a*a/5040;

When I've split the sum into steps, the model began to run. I guess the conclusion is that we can not make too many arithmetical operations for one "=" . The most interesting thing is that these lines stood after true if-statement and should not be even executed. Possible explanation is that some compilers start to process code after if-statement during its slow comparison. In our case it could start and then broke down on these long expressions.

I've implemented online frequency domain filter and applied it to MC_F.

+

Magnitude of the filter output at 1 Hz is the same as MC_F. This means that it is not hard for FIR to match the resonance. The problem is with the phase. We can not match the resonance exactly. If the resonance is at f₀ and we match at f₀ +/- df then in the frequency range (f₀, f₀ +/- df) the phase is not matched for 180. I guess the filter does not diverge because df is small but also the filter can not account for this huge phase lag. We need to slightly change the simulated actuator TF and see how the filter will react.

I'm trying to get some intuition how digital noise due to quantization shows up in iir filters. I decided to do tests in C using Python to calculate psd and visualize. I've implemented Direct Form 1, 2, "Biquad" and "Low Noise" forms of realization of second-order iir filter from Matt's presentation. There is a typo in the "Low Noise Form" scheme - a₁ and a₂ gains should be switched. Other then that schemes correctly implement 2 order iir. 

The input signal to each filter was a sine wave plus white noise with small amplitude x[n] = sin(2*pi*f*t[n]) + g*random( [-1, 1] ),  g << 1, f=1kHz. Sampling frequency was 16384 Hz. All 4 forms implemented 2 order low-pass butterworth filter with cut-off frequency 0.2 Hz

For g=1e-2 all implementations work fine. For g=1e-8 when quantization noise increases, all implementations give a lot of noise at low frequencies. I did not notice any significant difference between any of these implementations. I'll try to do more tests to figure out any difference in noise between the forms.

Quantization noise depends on the architecture of the processor, compiler and what not. But I do not think this can give a huge difference in results. We need to understand carefully digital noise during PSD estimation and all operations done at Matlab or Python.

 Direct Form 2 is noisy in the first test. This is the one similar to Matt's in his presentation. Input signal was a sine wave at 1 Hz with small amplitude white noise x[n] = sin(2*pi*1*t[n]) + 1e-10 * random( [-1, 1] ). It was filtered with a notch filter: f=1Hz, Q=100, depth=210dB. SOS representation was calculated in Foton. Sampling frequency is 16kHz.

DF2 output noise level is the same if I change white noise amplitude while DF1, BQF, LNF can follow it. Time series show quantization noise of DF2. I've plotted coherence of the signals relative to DF1, because non of the signals will be coherent to it at low frequencies due to fft calculations.  

In the second test the input was white noise  x[n] = random( [-1, 1] ) It was filtered with a 2 order low-pass butterworth filter with cut-off frequency f = 0.25 Hz. SOS representation was calculated in Python. Sampling frequency is 16kHz.

In this test all implementations work fine. I guess dtt works with single precision and for that reason we see disturbance in coherence when we do the same test online.

 This is because when we plot signals with sampling frequencies 2k and 16k with the same BW, we actually create psd/coherence using different numbers of points in FFT calculations as NFFT ~ fs/bw, fs-sampling frequency. So we secretly used 8 times more fft points for 16k signal then for 2k signal. Following plots demonstrate this effect. The first plot shows transfer function and coherence for filtering of 16k signal with butter('LowPass',4,8) and 2k signal with butter('LowPass',4,1)  when BW=0.1. There is a disturbance in coherence for 2k signal below 2 Hz. Now let's take BW=0.8 and plot transfer function and coherence for 16k signal. We can see the same effect of coherence disturbance.

The similar effect takes place when we change the cut-off frequency. The following plots show transfer function and coherence of two pairs of 2kHz signals. 4 order butterworth low-pass filter was used. For the first pair cut-off frequency was 1 Hz, for the second 10 Hz.  On the first plot BW=0.1 and there is a disturbance in coherence below 1 Hz. However on the second plot when BW=0.01, this effect is lost.

I guess my goal is to figure out when these effects come from fft calculations and when from digital filter noise.

 I wanted to test biquad form in this model. I added biquad=1 flag to cdsParameters, compiled, installed and restarted it. After that c1iscey suspended.

The same thing as we had several month ago

controls@c1iscey /opt/rtcds/caltech/c1/target/c1tst/c1tstepics 0$ cat iocC1.log

Starting iocInit
iocRun: All initialization complete
sh: iniChk.pl: command not found
Failed to load DAQ configuration file

 I've filtered a 1 Hz sin wave excitation with a notch filter inside c1sus and c1rfm models. The biquad key is switched on in the last one, c1sus uses DF2. The results are indeed different.

Still I do not like huge (2n+1) harmonics in the output of the biquad filter, I do not get them in the simulations. They are absent in the time series as well. So this is not a psd-estimation effect.

 Excitation generator created these harmonics. When I applied low-pass butterworth filter, I've got the result of online filtering close to simulations. On the second graph blue is biquad filter output spectrum, red corresponds to DF2. 1 Hz sin wave was filtered with a notch filter of Q=100, depth=300 at 1 Hz.

 We decided to write a script that will check online filters for digital noise. One method can be implemented using the following algorithm:

• calculate filter output using single precision
• calculate filter output using double precision and assume that it is precise
• find digital noise at the output of the filter when single precision is used
• extrapolate the result to the double precision filter dividing by 2^D-S ~ 10⁷, D - number of bits used in double precision mantissa, S - in single precision

Restriction: Single precision filter internal variables must be checked for overflows.

I applied this method to filtering a 1 Hz sine wave with a notch filter. Precise output should also be a 1 Hz wave => at other frequencies we see noise => digital noise spectrum should coincide with filter output. The plot shows the method worked out for this example.

Using this method I estimated digital noise of butter("LowPass", 2, 0.001) applied to white noise. Sampling frequency was 16 kHz. 

We have a Trillium for several days from Vladimir. I've put seismometer inside the foam box on linolium. I was not able to level the seismometer on granite as this Trillium does not have level screws. Does anybody know where they are?  Readout box stands on the foam box as seismometer cable is short (~2 meters).

Cables go to STS1 inputs (7-9) on ADC 3.

As we do not have legs for Trillium, I was advised to use shims to adjust the levels. However, they produce extra resonance at ~30 Hz + harmonics. Coherence is lost at these frequencies.

The script estimates digital noise produces by online filters. First version of Matlab files and complied c files are in scripts/digital_noise directory.

Algorithm for 1 filter bank (max number of filters = 10):

1.        extract sos - representation from Foton file for each filter (Matlab)
2.        download data from corresponding DQ channel using NDS (Matlab)
3.        find filters that are switched on (Matlab)
4.        filter signal using Df2 and BQF with single and double precision (C)
5.        estimate digital filter noise (Matlab)
6.        calculate power spectral density and plot the result (Matlab)

More details on (2)

Often DQ channels have reduced sampling rate. In this case the script will upsample data adding zeros.

AI filter is not applied. But in the end only the frequency range (0, DQ RATE / 2) is analyzed.

More details on (3):

• input switch ON / OFF => 4 / 0
• filters 1 - 6 ON /OFF
  • 1 => 48 / 0
  • 2 => 192 / 0
  • 3 => 768 / 0
  • 4 => 3072 / 0
  • 5 => 12288 / 0
  • 6 => 49152 / 0
  • If a switch is ON but there is no corresponding filter (one green and one red line under the switch) then the switch value is divided by 3

_SW2R channel value is the sum of the following numbers:

• decimation switch ON / OFF => 512 / 0
• output switch ON / OFF => 1024 / 0
• filters 7 - 10 ON /OFF
  
  • 7 => 3 / 0
  • 8 => 12 / 0
  • 9 => 48 / 0
  • 10 => 192 / 0
  • If a switch is ON but there is no corresponding filter (one green and one red line under the switch) then the switch value is divided by 3

Note: as for now Matlab script assumes that input, output and decimation filters are switched ON and there are no turned ON filter switches that do not correspond to any filters

More details on (5)

Digital noise using double precision is estimated by extrapolation of digital noise with single precision. The last is calculated by subtracting outputs of the filters with single and double precision. Then this noise is multiplied by 3 * 10^-7.

This extrapolation number was achieved by printf tests of the number 0.123456789012345678 with single and double precision on C. Using type 'float' variables 10 significant numbers show up, using type 'double' - 17.

I also did 'calibration tests' to achieve extrapolation number - signal was filters with an aggresive low-pass filter. At high frequencies filter output spectrum is flat => digital noise amplitude must be the same. The plot shows GUR1_X channel filtered with low-pass chebyshev type 1 filter.

However, extrapolation number is not the same for all cases. In the following example of analyzing BS_SUSPOS filter bank using extrapolation 3 * 10^-7 we get noise that is slightly overestimated. In some other examples we need to take a larger number. But in average, I think, this is a good approximation.

To avoid extrapolation problem we can use long double precision (~19 digits). I was able to do this with gcc compiler. However, in mex compiler using long double in filter calculations, I do not get any better precision then using double precision. I'll think more about it.

 After that RFM -> OAF communication through PCIE became bad again. Inside CommData2.c cache flushing is not allowed

// If PCIE comms show errors, may want to add this cache flushing
            #if 0
            if(ipcInfo[ii].netType == IPCIE)
                clflush_cache_range (ipcInfo[ii].pIpcData->dBlock[sendBlock][ipcIndex].data, 16);
            #endif

As a result, a significant part of MC_F and other signals is lost during RFM -> OAF transmission (270 - 330 out of 2048 per second)

Last time when I replaced 0 for 1, it suspended SUS machine because of the code bug. Alex modified a couple of files in the old version and it started to work. Do you know if this bug is fixed in the new version?

 All optics except MC2 and ETMX are crazy

Alex also modified RCG script to generate -O in the Makefile for c1pem model:

controls@pianosa:/opt/rtcds/rtscore/release/src/epics/util 127$ svn diff
feCodeGen.pl 
Index: feCodeGen.pl
===================================================================
--- feCodeGen.pl (revision 2999)
+++ feCodeGen.pl (working copy)
@@ -3183,7 +3183,12 @@
print OUTM "\n";
}
print OUTM "ALL \+= user_mmap \$(TARGET_RTL)\n";
+# do not optimize c1pem
+if ($skeleton eq "c1pem") {
+print OUTM "EXTRA_CFLAGS += -O -w -I../../include\n";
+} else {
print OUTM "EXTRA_CFLAGS += -O3 -w -I../../include\n";
+}
print OUTM "EXTRA_CFLAGS += -I/opt/gm/include\n";
print OUTM "EXTRA_CFLAGS += -I/opt/mx/include\n";

 Using Guralp, STS-2 and Trillium I compared Gur and STS-2 self-noise assuming that Trillium noise is not worse then STS-2 noise.

Interesting that STS-2 (or Trillium if its noise is worse) noise is not too much better then Guralp noise.

STS-2 - end of X arm

GUR 2 - isolation box

TRILLIUM - 1Y3 (DC power supply uses 1Y3 AC power, please do not close the door completely)

GUR 1 - end of Y arm

Now we have several "triangular seismic antennas". Different configurations can be chosen to compare the results.

 Probably, security. You can call 5555 and ask them. Otherwise you can ask them to come and check everything.

 A huge data flow goes from PEM to OAF through RFM. I think we need to make PEM and OAF run on the same machine and transmit signals through the shared memory.

A Matlab script to calculate Wiener filter coefficients and convert fir to iir is ready. Input is a file with zero mean witness and desired signals, output is a Foton zpk command to specify iir filter.

The plot shows comparison of offline fir , iir and online iir filtering. Spectrum below 4 Hz is still oscillating due to acoustic coupling, this is not a filtering effect. At 1 Hz actuator is badly compensated, more work should be done. Other then that online and offline filtering are the same.

The lack of AA filter for MCL signal is RFM model strongly disturbed entering to OAF signal

Transmission of signals between RFM and OAF is bad again. Now we do not see any errors in IPC_ERR monitors so models think that they get all data but the data is wrong

I've put EM 172 microphones inside Steve's isolation box to measure their noise. I've attached mics to each other and aligned them using the tape.

At low frequencies (below 1 Hz) the noise is limited by ADC as there is a 10 Hz high-pass filter inside mic readout box.

ADC noise is measured by splitting the signal from 1 mic into 2 ADC channels.

Data from PEM now goes directly to OAF without using RFM. Transmission RFM -> OAF errors are gone as RFM has to read 30 channels less now.

Again kernel "protection error" occured as before with PEM model so OAF model could not start. I changed optimization flag to -02, this fixed the problem.

Static filter was adjusted to filter 1 Hz resonance in MCL and it could do it. Stack is not great in this experiment due to the phase mismatch. I'll fix it.

Jenne and I did adaptive filtering of MC_L and measured how X and Y ARM control signals change compared to non-filtered MC_L. We did the test during 1.5 Hz seismic noise activity and adaptive filter was able to subtract it. However, it adds noise at high frequencies, It is not seen in MC_L but it is present in the ARMs control signals.

I'll investigate this problem. May be we need to reduce adaptation gain. In this experiment it was 0.1 and adaptive filter convergence time was equal to 1-2 mins.

I did a simulation of linear quadratic gaussian (LQG) controller applied to local damping. The cost function was frequency shaped to have a peak at 1 Hz. This technique prevents the controller from adding sensor noise at high and very low frequencies.

Noise was simulated to have 1/f spectrum (seismic) multiplied by stack with a resonance at 4 Hz with Q=5.

I did a simulation of a cavity, feedback signal was calculated using LQG controller. I assumed that there is not length -> angle coupling and 2 mirrors that form the cavity have the same equation of motions (Q and eigen frequencies are the same). Cost functional was chosen in such a way that frequencies below 15 Hz contribute much more then other frequencies.

Gains in the controller are calculated to minimize the cost functional.

This technique works well, but it requires full information about the system states. If we do not assume that cavity mirrors have the same equations of motion then we need to apply Kalman filter to approximate the position of one of the mirrors.

 I've installed the mirrors on 4 tip-tilts. I was able to align 3 of them in pitch, the last one has a screw with damaged thread, I'll continue with it tomorrow.

Alignment accuracy in pitch is ~0.1 mrad. Mirrors oscillate a lot probably due to air flow coming from the side wall.

I've measured seismic and acoustic noise on BS and AS tables. It seems that horizontal motion of BS table is ~1.5-2 times more then AS table in the frequency range 5-50 Hz.

Edit by Den: this was POI table, not BS!

I've suspended microphones around the lab

 This seismic measurement is for BS and AS tables.

Time domain control using LQR technique is now applied to ETMX sus position. The plan was to do it for oplevs, I'll do it after the vent.

The cost function for state space variables was determined by TF  900 / (s + 30)^2. There was no penulty imposed for velocity, only for position. We can try that configuration as well.

 We would appreciate some plots. Learning curves of recurrent NN working as a plant are interesting. For harmonic oscillator your RNN should not contain any hidden layers - only 1 input and 1 output node and 2 delays at each of them. Activation function should be linear. If your code is correct, this configuration will match oscillator perfectly. The question is how much time does it take to adapt.

Does FANN support regularization? I think this will make your controller more stable. Try to use more advanced algorithms then gradient descent for adaptation. They will increase convergence speed. For example, look at fminunc function at Matlab.

I did TT alignment using red laser and QPD.

 I had a problem aligning TT with frame number SN-035 as some screws are damaged so all what I could move were 2 blades on the sides of the mirrors. But this was not enough to align pitch and yaw simultaneously.  It is possible to align pitch only, but then I got a huge yaw angle (~0.05-0.1 rad). The only option I had was to make a reasonable alignment in yaw and then suspend several washes on the screw on the bottom of the mirror to align pitch.

Attached are flag positions inside coils. 1 - SN034, 2 - SN012, 3 - SN006, 4 - SN035. For each TT there are 4 pictures with flag 1-4: UL, LL, UR, LR

I've installed Guralp readout box back and it turned out that it does not work with voltage provided from the rack (+13.76 0 -14.94).  +/-12 voltage regulators inside the box convert it to -0.9 0 -12. I've connected the box to +/-15 DC voltage supply to measure seismic motion at the ETMY table. Readout box works fine with +/- 15.

Seismic noise on the ETMY table measured to be a few times higher then on the floor in horizontal direction in the frequency range 50 - 200 Hz. Attached are compared spectrums of X, Y and Z motions.

I've applied LQR approach to MC_L locking. Results show that LQR does not make MC_F signal smaller below 0.3 Hz in contrast with classical locking. This might indicate that in this frequency range we see sensing noise as LQR was provided with state-space model of MC only so it tries to reduce displacement noise. It is also possible that state-space model is not accurate enough.

  [Koji, Steve, Den]

TT alignment is fine, yaw damping is satisfactory, pitch damping is slow. We might want to add magnets to the mirror and attach blades to the frame for pitch edge current damping.

We are moving towards electronics testing.

PZT monitor is not lying to us, I've measured it with a voltmeter. But PMC SERVO is still interesting. If I break the loop after PMCERR signal monitor (C1:PSL-PMC_BLANK=0), I can change PZT voltage from 0V (C1:PSL-PMC_RAMP = -7.3) to 132V (C1:PSL-PMC_RAMP=-10.0). If I break the loop by enabling TEST1, PZT voltage goes up to 294 V, though voltage on the TEST1 and MIXER MON is 0.

 Do you know how to burtrestore or restart slow machines?

Edit by Den: I did burtrestore of c1psl.snap from 2 days ago. Still slow machines behave not normal. For example, if I sweep C1:PSL-FSS_SLOWDC, SLOW monitor value does not change.

Guralp readout box received +13.7 /0/ -15V instead +15V because of the broken fuse. Power provided by the source is normal.

Edit by Den: I've found a similar fuse on one of the tables and borrowed it. Guralp is not working again.

 I've added calibration gains to Guralp (to um/sec) and EM172 (to Pa) channels.

We can run PEM at 16 kHz. I think Foton file stores both sos-representation and filter commands which are independent of the sampling frequency, so it should be possible to change model sampling frequency quickly.

In fact, we can save data at 64 kHz from iop models. I've done this once with MC_F channel. However, I did not test EM172 noise at frequencies > 1 kHz.

 I've meant Guralp is NOW working again 

 MC2 was tripped again. I think the answer is that watchdog's critical value was too small C1:SUS-MC2_PD_MAX_VAR = 10, so seismic could trip MC2. I've changed the value to 100.

I've applied online state estimation technique using Kalman filter to LQG controller. It helps to estimate states that we do not measure. I've considered MC2 local damping, we measure position and want to estimate velocity that we need for control. We can either differentiate the signal or apply state estimation to avoid huge noise injection at high frequencies. In state estimation we need to know noise covariance, I've assumed that LID sensor noise is 0.1 nm. Though covariance can be calculated better.

In the time-domain figure C1:SUS-MC2_SUSPOS_IN1 = MC2 postion, C1:SUS-MC2_SUSPOS_OUT = MC2 velocity obtained by differentiation, 2 other channels are estimations of position and velocity.

 I think state estimation and optimal control are two different techniques that are often used together. Sometimes (as for pendulum) we can use LQG without state estimation as we need only position and velocity. But for more complex systems (like quad suspension) the states of all 4 masses can be reconstructed in some optimal way using information from only one of them if the dynamics is sufficiently well known. When current system states are measured/estimated we can apply control where all our filters are hidden.

 The algorithm needs to know about expected seismic noise transmission from the ground to the optic, but it might be not very precise. I gave it some rough estimate, there are better ways to do it. I think that we'll understand whether we need state estimation or not when we'll move to more complex systems. Brett uses a similar approach for his modal control. Interesting if these methods + seismometer readings will be able to say if one of your sensors is noisier then others.

I went inside to align the beam on WFS and noticed that oscillations in yaw are ~10 times stronger then in pitch. I've plot rms of pitch and yaw measured by LID sensors and saw that MC3 yaw rms motion is a few times larger then pitch.

Also MC1 input diag matrix does not diagonalize signals for pitch and yaw. In the spectrums of these signals all 4 resonance are equally seen during the free swinging. I think we should rediagonalize MC1.

Another thing is that if MC1 and MC3 are on the same stack,  pitch and yaw spectrums of these mirrors should be comparable. But MC1 signal is ~2-3 times larger then of MC3. I think we should correct calibration.