40m QIL Cryo_Lab CTN SUS_Lab TCS_Lab OMC_Lab CRIME_Lab FEA ENG_Labs OptContFac Mariner WBEEShop
  40m Log, Page 296 of 330  Not logged in ELOG logo
ID Date Authordown Type Category Subject
  16151   Fri May 21 09:44:52 2021 Ian MacMillanUpdateCDSSUS simPlant model

The transfer function given in the previous post was slightly incorrect the units did not make sense the new function is:

\frac{x}{F}=\frac{1}{m\omega_0^2-m\omega^2+im\frac{\omega_0 \omega }{Q}}

I have attached a quick derivation below in attachment 1

Attachment 1: Transfer_Function_of_Damped_Harmonic_Oscillator.pdf
Transfer_Function_of_Damped_Harmonic_Oscillator.pdf
  16153   Fri May 21 14:36:20 2021 Ian MacMillanUpdateCDSSUS simPlant model

The plant transfer function of the pendulum in the s domain is:

H(s)=\frac{x(s)}{F(s)}=\frac{1}{ms^2+m\frac{\omega_0}{Q}s+m\omega_0^2}

Using Foton to make a plot of the TF needed and using m=40kg, w0=3Hz, and Q=50 (See attachment 1). It is easiest to enter the above filter using RPoly and saved it as Plant_V1

Attachment 1: Plant_Mod_TF.pdf
Plant_Mod_TF.pdf
  16177   Thu Jun 3 13:06:47 2021 Ian MacMillanUpdateCDSSUS simPlant model

I was able to measure the transfer function of the plant filter module from the channel X1:SUP-C1_SUS_SINGLE_PLANT_Plant_POS_Mod_EXC to X1:SUP-C1_SUS_SINGLE_PLANT_Plant_POS_Mod_OUT. The resulting transfer function is shown below. I have also attached the raw data for making the graph.

Next, I will make a script that will make the photon filters for all the degrees of freedom and start working on the matrix version of the filter module so that there can be multiple degrees of freedom.

Attachment 1: SingleSusPlantTF.pdf
SingleSusPlantTF.pdf
Attachment 2: SUS_PLANT_TF.zip
  16191   Mon Jun 7 17:49:19 2021 Ian MacMillanUpdateCDSSUS simPlant model

Added difference to the graph. I included the code so that others could see what it looks like and use it for easy use.

Attachment 1: SingleSusPlantTF.pdf
SingleSusPlantTF.pdf
Attachment 2: TF_Graph_Code.zip
  16195   Wed Jun 9 13:50:48 2021 Ian MacMillanUpdateCDSSUS simPlant model

I have attached an updated transfer function graph with the residual easier to see. I thought here I would include a better explanation of what this transfer function was measuring.

This transfer function was mainly about learning how to use DTT and Foton to make and measure transfer functions. Therefore it is just measuring across a single CDS filter block. X1SUP_ISOLATED_C1_SUS_SINGLE_PLANT_Plant_POS_Mod block to be specific. This measurement shows that the block is doing what I programmed it to do with Foton. The residual is probably just because the measured TF had fewer points than the calculated one.

The next step is to take a closed-loop TF of the system and the control module.

After that, I want to add more degrees of freedom to the model. both in the plant and in the controls.

Attachment 1: SingleSusPlantTF.pdf
SingleSusPlantTF.pdf
  16201   Tue Jun 15 11:46:40 2021 Ian MacMillanUpdateCDSSUS simPlant model

I have added more degrees of freedom. The model includes x, y, z, pitch, yaw, roll and is controlled by a matrix of transfer functions (See Attachment 2). I have added 5 control filters to individually control UL, UR, LL, LR, and side. Eventually, this should become a matrix too but for the moment this is fine.

Note the Unit delay blocks in the control in Attachment 1. The model will not compile without these blocks.

Attachment 1: x1sup_isolated-6-15-v1.pdf
x1sup_isolated-6-15-v1.pdf
Attachment 2: C1_SUS_SINGLE_PLANT-6-15-v1.pdf
C1_SUS_SINGLE_PLANT-6-15-v1.pdf
  16208   Thu Jun 17 11:19:37 2021 Ian MacMillanUpdateCDSCDS Upgrade

Jon and I tested the ADC and DAC cards in both of the systems on the test stand. We had to swap out an 18-bit DAC for a 16-bit one that worked but now both machines have at least one working ADC and DAC.

[Still working on this post. I need to look at what is in the machines to say everything ]

  16217   Mon Jun 21 17:15:49 2021 Ian MacMillanUpdateCDSCDS Upgrade

Anchal and I wrote a script (Attachment 1) that will test the ADC and DAC connections with inputs on the INMON from -3000 to 3000. We could not run it because some of the channels seemed to be frozen. 

Attachment 1: DAC2ADC_Test.py
import os
import time
import numpy as np
import subprocess
from traceback import print_exc
import argparse


def grabInputArgs():
    parser = argparse.ArgumentParser(
... 75 more lines ...
  16220   Tue Jun 22 16:53:01 2021 Ian MacMillanUpdateCDSFront-End Assembly and Testing

The channels on both the C1BHD and C1SUS2 seem to be frozen: they arent updating and are holding one value. To fix this Anchal and I tried:

  • restarting the computers 
    • restarting basically everything including the models
  • Changing the matrix values
  • adding filters
  • messing with the offset 
  • restarting the network ports (Paco suggested this apparently it worked for him at some point)
  • Checking to make sure everything was still connected inside the case (DAC, ADC, etc..)

I wonder if Jon has any ideas. 

  16224   Thu Jun 24 17:32:52 2021 Ian MacMillanUpdateCDSFront-End Assembly and Testing

Anchal and I ran tests on the two systems (C1-SUS2 and C1-BHD). Attached are the results and the code and data to recreate them.

We connected one DAC channel to one ADC channel and thus all of the results represent a DAC/ADC pair. We then set the offset to different values from -3000 to 3000 and recorded the measured signal. I then plotted the response curve of every DAC/ADC pair so each was tested at least once.

There are two types of plots included in the attachments

1) a summary plot found on the last pages of the pdf files. This is a quick and dirty way to see if all of the channels are working. It is NOT a replacement for the other plots. It shows all the data quickly but sacrifices precision.

2) In an in-depth look at an ADC/DAC pair. Here I show the measured value for a defined DC offset. The Gain of the system should be 0.5 (put in an offset of 100 and measure 50). I included a line to show where this should be. I also plotted the difference between the 0.5 gain line and the measured data. 

As seen in the provided plots the channels get saturated after about the -2000 to 2000 mark, which is why the difference graph is only concentrated on -2000 to 2000 range. 

Summary: all the channels look to be working they all report very little deviation off of the theoretical gain. 

Note: ADC channel 31 is the timing signal so it is the only channel that is wildly off. It is not a measurement channel and we just measured it by mistake.

Attachment 1: C1-SU2_Channel_Responses.pdf
C1-SU2_Channel_Responses.pdf C1-SU2_Channel_Responses.pdf C1-SU2_Channel_Responses.pdf C1-SU2_Channel_Responses.pdf C1-SU2_Channel_Responses.pdf C1-SU2_Channel_Responses.pdf C1-SU2_Channel_Responses.pdf C1-SU2_Channel_Responses.pdf
Attachment 2: C1-BHD_Channel_Responses.pdf
C1-BHD_Channel_Responses.pdf C1-BHD_Channel_Responses.pdf C1-BHD_Channel_Responses.pdf C1-BHD_Channel_Responses.pdf C1-BHD_Channel_Responses.pdf C1-BHD_Channel_Responses.pdf C1-BHD_Channel_Responses.pdf C1-BHD_Channel_Responses.pdf
Attachment 3: CDS_Channel_Test.zip
  16230   Wed Jun 30 14:09:26 2021 Ian MacMillanUpdateCDSSUS simPlant model

I have looked at my code from the previous plot of the transfer function and realized that there is a slight error that must be fixed before we can analyze the difference between the theoretical transfer function and the measured transfer function.

The theoretical transfer function, which was generated from Photon has approximately 1000 data points while the measured one has about 120. There are no points between the two datasets that have the same frequency values, so they are not directly comparable. In order to compare them I must infer the data between the points. In the previous post [16195] I expanded the measured dataset. In other words: I filled in the space between points linearly so that I could compare the two data sets. Using this code:

#make values for the comparison
tck_mag = splrep(tst_f, tst_mag) # get bspline representation given (x,y) values
gen_mag = splev(sim_f, tck_mag) # generate intermediate values
dif_mag=[]
for x in range(len(gen_mag)):
    dif_mag.append(gen_mag[x]-sim_mag[x]) # measured minus predicted

tck_ph = splrep(tst_f, tst_ph) # get bspline representation given (x,y) values
gen_ph = splev(sim_f, tck_ph) # generate intermediate values
dif_ph=[]
for x in range(len(gen_ph)):
    dif_ph.append(gen_ph[x]-sim_ph[x])

At points like a sharp peak where the measured data set was sparse compared to the peak, the difference would see the difference between the intermediate “measured” values and the theoretical ones, which would make the difference much higher than it really was.

To fix this I changed the code to generate the intermediate values for the theoretical data set. Using the code here:

tck_mag = splrep(sim_f, sim_mag) # get bspline representation given (x,y) values
gen_mag = splev(tst_f, tck_mag) # generate intermediate values
dif_mag=[]
for x in range(len(tst_mag)):
    dif_mag.append(tst_mag[x]-gen_mag[x])#measured minus predicted

tck_ph = splrep(sim_f, sim_ph) # get bspline representation given (x,y) values
gen_ph = splev(tst_f, tck_ph) # generate intermediate values
dif_ph=[]
for x in range(len(tst_ph)):
    dif_ph.append(tst_ph[x]-gen_ph[x])

Because this dataset has far more values (about 10 times more) the previous problem is not such an issue. In addition, there is never an inferred measured value used. That makes it more representative of the true accuracy of the real transfer function.

This is an update to a previous plot, so I am still using the same data just changing the way it is coded. This plot/data does not have a Q of 1000. That plot will be in a later post along with the error estimation that we talked about in this week’s meeting.

The new plot is shown below in attachment 1. Data and code are contained in attachment 2

Attachment 1: SingleSusPlantTF.pdf
SingleSusPlantTF.pdf
Attachment 2: Plant_TF_Test.zip
  16289   Mon Aug 23 15:25:59 2021 Ian MacMillanUpdateCDSSUS simPlant model

I am adding a State-space block to the SimPlant cds model using the example Chris gave. I made a new folder in controls called SimPlantStateSpace. wI used the code below to make a state-space LTI model with a 1D pendulum then I converted it to a discrete system using c2d matlab function. Then I used these in the rtss.m file to create the state space code I need in the SimPlantStateSpace_1D_model.h file.

%sys_model.m

Q = 1000;
phi = 1/Q;
g = 9.806;
m = 0.24; % mass of pendulum
l = 0.248; %length of pendulum
w_0 = sqrt(g/l);

f=16000 %this is the frequency of the channel that will be used

A = [0 1; -w_0^2*(1+1/Q*1i) -w_0/Q]
B = [0; 1/m];
C = [1 0];
D = [0];
sys_dc = ss(A,B,C,D)

sys=c2d(sys_dc, 1/f)

This code outputs the discrete state space that is added to the header file attached.

Attachment 1: SimPlantStateSpace.zip
  16339   Thu Sep 16 14:08:14 2021 Ian MacMillanFrogs Tour

I gave some of the data analysts a look around because they asked and nothing was currently going on in the 40m. Nothing was changed.

  16348   Mon Sep 20 15:42:44 2021 Ian MacMillanSummaryComputersQuantization Code Summary

This post serves as a summary and description of code to run to test the impacts of quantization noise on a state-space implementation of the suspension model.

Purpose: We want to use a state-space model in our suspension plant code. Before we can do this we want to test to see if the state-space model is prone to problems with quantization noise. We will compare two models one for a standard direct-ii filter and one with a state-space model and then compare the noise from both. 

Signal Generation:

First I built a basic signal generator that can produce a sine wave for a specified amount of time then can produce a zero signal for a specified amount of time. This will let the model ring up with the sine wave then decay away with the zero signal. This input signal is generated at a sample rate of 2^16 samples per second then stored in a numpy array. I later feed this back into both models and record their results.

State-space Model:

The code can be seen here

The state-space model takes in the list of excitation values and feeds them through a loop that calculates the next value in the output.

Given that the state-space model follows the form

  \dot{x}(t)=\textbf{A}x(t)+ \textbf{B}u(t)   and  y(t)=\textbf{C}x(t)+ \textbf{D}u(t) ,

the model has three parts the first equation, an integration, and the second equation.

  1. The first equation takes the input x and the excitation u and generates the x dot vector shown on the left-hand side of the first state-space equation.
  2. The second part must integrate x to obtain the x that is found in the next equation. This uses the velocity and acceleration to integrate to the next x that will be plugged into the second equation
  3. The second equation in the state space representation takes the x vector we just calculated and then multiplies it with the sensing matrix C. we don't have a D matrix so this gives us the next output in our system

This system is the coded form of the block diagram of the state space representation shown in attachment 1

Direct-II Model:

The direct form 2 filter works in a much simpler way. because it involves no integration and follows the block diagram shown in Attachment 2, we can use a single difference equation to find the next output. However, the only complication that comes into play is that we also have to keep track of the w(n) seen in the middle of the block diagram. We use these two equations to calculate the output value

y[n]=b_0 \omega [n]+b_1 \omega [n-1] +b_2 \omega [n-2],  where w[n] is  \omega[n]=x[n] - a_1 \omega [n-1] -a_2 \omega[n-2]

Bit length Control:

To control the bit length of each of the models I typecast all the inputs using np.float then the bit length that I want. This simulates the computer using only a specified bit length. I have to go through the code and force the code to use 128 bit by default. Currently, the default is 64 bit which so at the moment I am limited to 64 bit for highest bit length. I also need to go in and examine how numpy truncates floats to make sure it isn't doing anything unexpected.

Bode Plot: 

The bode plot at the bottom shows the transfer function for both the IIR model and the state-space model. I generated about 100 seconds of white noise then computed the transfer function as 

G(f) = \frac{P_{csd}(f)}{P_{psd}(f)}

which is the cross-spectral density divided by the power spectral density. We can see that they match pretty closely at 64 bits. The IIR direct II model seems to have more noise on the surface but we are going to examine that in the next elog

 

Attachment 1: 472px-Typical_State_Space_model.svg.png
472px-Typical_State_Space_model.svg.png
Attachment 2: Biquad_filter_DF-IIx.svg.png
Biquad_filter_DF-IIx.svg.png
Attachment 3: SS-IIR-TF.pdf
SS-IIR-TF.pdf
  16355   Wed Sep 22 14:22:35 2021 Ian MacMillanSummaryComputersQuantization Noise Calculation Summary

Now that we have a model of how the SS and IIR filters work we can get to the problem of how to measure the quantization noise in each of the systems. Den Martynov's thesis talks a little about this. from my understanding: He measured quantization noise by having two filters using two types of variables with different numbers of bits. He had one filter with many more bits than the second one. He fed the same input signal to both filters then recorded their outputs x_1 and x_2, where x_2 had the higher number of bits. He then took the difference x_1-x_2. Since the CDS system uses double format, he assumes that quantization noise scales with mantissa length. He can therefore extrapolate the quantization noise for any mantissa length.

Here is the Code that follows the following procedure (as of today at least)

This problem is a little harder than I had originally thought. I took Rana's advice and asked Aaron about how he had tackled a similar problem. We came up with a procedure explained below (though any mistakes are my own):

  1. Feed different white noise data into three of the same filter this should yield the following equation: \textbf{S}_i^2 =\textbf{S}_{ni}^2+\textbf{S}_x^2, where \textbf{S}_i^2 is the power spectrum of the output for the ith filter,  \textbf{S}_{ni}^2 is the noise filtered through an "ideal" filter with no quantization noise, and  \textbf{S}_x^2 is the power spectrum of the quantization noise. Since we are feeding random noise into the input the power of the quantization noise should be the same for all three of our runs.
  2. Next, we have our three outputs:  \textbf{S}_1^2,  \textbf{S}_2^2, and  \textbf{S}_3^2 that follow the equations: 

\textbf{S}_1^2 =\textbf{S}_{n1}^2+\textbf{S}_x^2

\textbf{S}_2^2 =\textbf{S}_{n2}^2+\textbf{S}_x^2

\textbf{S}_3^2 =\textbf{S}_{n3}^2+\textbf{S}_x^2

From these three equations, we calculate the three quantities: \textbf{S}_{12}^2\textbf{S}_{23}^2, and \textbf{S}_{13}^2 which are calculated by:

\textbf{S}_{12}^2 =\textbf{S}_{1}^2-\textbf{S}_2^2\approx \textbf{S}_{n1}^2 -\textbf{S}_{n2}^2

\textbf{S}_{23}^2 =\textbf{S}_{2}^2-\textbf{S}_3^2\approx \textbf{S}_{n2}^2 -\textbf{S}_{n3}^2

\textbf{S}_{13}^2 =\textbf{S}_{1}^2-\textbf{S}_3^2\approx \textbf{S}_{n1}^2 -\textbf{S}_{n3}^2

from these quantities, we can calculate three values: \bar{\textbf{S}}_{n1}^2\bar{\textbf{S}}_{n2}^2, and \bar{\textbf{S}}_{n3}^2 since these are just estimates we are using a bar on top. These are calculated using:

\bar{\textbf{S}}_{n1}^2\approx\frac{1}{2}(\textbf{S}_{12}^2+\textbf{S}_{13}^2+\textbf{S}_{23}^2)

\bar{\textbf{S}}_{n2}^2\approx\frac{1}{2}(-\textbf{S}_{12}^2+\textbf{S}_{13}^2+\textbf{S}_{23}^2)

\bar{\textbf{S}}_{n3}^2\approx\frac{1}{2}(\textbf{S}_{12}^2+\textbf{S}_{13}^2-\textbf{S}_{23}^2)

using these estimates we can then estimate  \textbf{S}_{x}^2  using the formula:

\textbf{S}_{x}^2 \approx \textbf{S}_{1}^2 - \bar{\textbf{S}}_{n1}^2 \approx \textbf{S}_{2}^2 - \bar{\textbf{S}}_{2}^2 \approx \textbf{S}_{3}^2 - \bar{\textbf{S}}_{n3}^2

we can average the three estimates for  \textbf{S}_{x}^2  to come up with one estimate.

This procedure should be able to give us a good estimate of the quantization noise. However, in the graph shown in the attachments below show that the noise follows the transfer function of the model to begin with. I would not expect this to be true so I believe that there is an error in the above procedure or in my code that I am working on finding. I may have to rework this three-corner hat approach. I may have a mistake in my code that I will have to go through.

I would expect the quantization noise to be flatter and not follow the shape of the transfer function of the model. Instead, we have what looks like just the result of random noise being filtered through the model. 

Next steps:

The first real step is being able to quantify the quantization noise but after I fix the issues in my code I will be able to start liking at optimal model design for both the state-space model and the direct form II model. I have been looking through the book "Quantization noise" by Bernard Widrow and Istvan Kollar which offers some good insights on how to minimize quantization noise. 

Attachment 1: IIR64-bitnoisespectrum.pdf
IIR64-bitnoisespectrum.pdf
  16360   Mon Sep 27 12:12:15 2021 Ian MacMillanSummaryComputersQuantization Noise Calculation Summary

I have not been able to figure out a way to make the system that Aaron and I talked about. I'm not even sure it is possible to pull the information out of the information I have in this way. Even the book uses a comparison to a high precision filter as a way to calculate the quantization noise:

"Quantization noise in digital filters can be studied in simulation by comparing the behavior of the actual quantized digital filter with that of a refrence digital filter having the same structure but whose numerical calculations are done extremely accurately."
-Quantization Noise by Bernard Widrow and Istvan Kollar (pg. 416)

Thus I will use a technique closer to that used in Den Martynov's thesis (see appendix B starting on page 171). A summary of my understanding of his method is given here:

A filter is given raw unfiltered gaussian data f(t) then it is filtered and the result is the filtered data x(t) thus we get the result: f(t)\rightarrow x(t)=x_N(t)+x_q(t)  where x_N(t) is the raw noise filtered through an ideal filter and x_q(t) is the difference which in this case is the quantization noise. Thus I will input about 100-1000 seconds of the same white noise into a 32-bit and a 64-bit filter. (hopefully, I can increase the more precise one to 128 bit in the future) then I record their outputs and subtract the from each other. this should give us the Quantization error e(t):
e(t)=x_{32}(t)-x_{64}(t)=x_{N_{32}}(t)+x_{q_{32}}(t) - x_{N_{64}}(t)-x_{q_{64}}(t)
and since x_{N_{32}}(t)=x_{N_{64}}(t) because they are both running through ideal filters:
e(t)=x_{N}(t)+x_{q_{32}}(t) - x_{N}(t)-x_{q_{64}}(t)
e(t)=x_{q_{32}}(t) -x_{q_{64}}(t)
and since in this case, we are assuming that the higher bit-rate process is essentially noiseless we get the Quantization noise x_{q_{32}}(t).

If we make some assumptions, then we can actually calculate a more precise version of the quantization noise:

"Since aLIGO CDS system uses double precision format, quantization noise is extrapolated assuming that it scales with mantissa length"
-Denis Martynov's Thesis (pg. 173)

From this assumption, we can say that the noise difference between the 32-bit and 64-bit filter outputs:  x_{q_{32}}(t)-x_{q_{64}}(t)  is proportional to the difference between their mantissa length. by averaging over many different bit lengths, we can estimate a better quantization noise number.

I am building the code to do this in this file

  16361   Mon Sep 27 16:03:15 2021 Ian MacMillanSummaryComputersQuantization Noise Calculation Summary

I have coded up the procedure in the previous post: The result does not look like what I would expect. 

As shown in Attachment1 I have the power spectrum of the 32-bit output and the 64-bit output as well as the power spectrum of the two subtracted time series as well as the subtracted power spectra of both. unfortunately, all of them follow the same general shape of the raw output of the filter. 

I would not expect quantization noise to follow the shape of the filter. I would instead expect it to be more uniform. If anything I would expect the quantization noise to increase with frequency. If a high-frequency signal is run through a filter that has high quantization noise then it will severely degrade: i.e. higher quantization noise. 

This is one reason why I am confused by what I am seeing here. In all cases including feeding the same and different white noise into both filters, I have found that the calculated quantization noise is proportional to the response of the filter. this seems wrong to me so I will continue to play around with it to see if I can gain any intuition about what might be happening.

Attachment 1: DeltaNoiseSpectrum.pdf
DeltaNoiseSpectrum.pdf
  16366   Thu Sep 30 11:46:33 2021 Ian MacMillanSummaryComputersQuantization Noise Calculation Summary

First and foremost I have the updated bode plot with the mode moved to 10 Hz. See Attachment 1. Note that the comparison measurement is a % difference whereas in the previous bode plot it was just the difference. I also wrapped the phase so that jumps from -180 to 180 are moved down. This eliminates massive jumps in the % difference. 

Next, I have two comparison plots: 32 bit and 64bit. As mentioned above I moved the mode to 10 Hz and just excited both systems at 3.4283Hz with an amplitude of 1. As we can see on the plot the two models are practically the same when using 64bits. With the 32bit system, we can see that the noise in the IIR filter is much greater than in the State-space model at frequencies greater than our mode.

Note about windowing and averaging: I used a Hanning window with averaging over 4 neighbor points. I came to this number after looking at the results with less averaging and more averaging. In the code, this can be seen as nperseg=num_samples/4 which is then fed into signal.welch

Attachment 1: SS-IIR-Bode.pdf
SS-IIR-Bode.pdf
Attachment 2: PSD_32bit.pdf
PSD_32bit.pdf
Attachment 3: PSD_64bit.pdf
PSD_64bit.pdf
  16403   Thu Oct 14 16:38:26 2021 Ian MacMillanUpdateGeneralKicking optics in freeSwing measurment

[Ian, Anchal]

We are going to kick the optics tonight at 2am.

The optics we will kick are the PRM BS ITMX ITMY ETMX ETMY

We will kick each one once and record for 2000 seconds and the log files will be placed in users/ian/20211015_FreeSwingTest/logs.

  16406   Fri Oct 15 12:14:27 2021 Ian MacMillanUpdateGeneralKicking optics in freeSwing measurment

[Ian, Anchal]

we ran the free swinging test last night and the results match up with in 1/10th of a Hz. We calculated the peak using the getPeakFreqs2 script to find the peaks and they are close to previous values from 2016.

In attachment 1 you will see the results of the test for each optic.

The peak values are as follows:

Optic POS (Hz) PIT (Hz) YAW (Hz) SIDE (Hz)
PRM 0.94 0.96 0.99 0.99
MC2 0.97 0.75 0.82 0.99
ETMY 0.98 0.98 0.95 0.95
MC1 0.97 0.68 0.80 1.00
ITMX 0.95 0.68 0.68 0.98
ETMX 0.96 0.73 0.85 1.00
BS 0.99 0.74 0.80 0.96
ITMY 0.98 0.72 0.72 0.98
MC3 0.98 0.77 0.84 0.97

The results from 2016 can be found at: /cvs/cds/rtcdt/caltech/c1/scripts/SUS/PeakFit/parameters2.m

Attachment 1: 20211015_Kicktest_plot.pdf
20211015_Kicktest_plot.pdf 20211015_Kicktest_plot.pdf 20211015_Kicktest_plot.pdf 20211015_Kicktest_plot.pdf 20211015_Kicktest_plot.pdf 20211015_Kicktest_plot.pdf 20211015_Kicktest_plot.pdf 20211015_Kicktest_plot.pdf
  16414   Tue Oct 19 18:20:33 2021 Ian MacMillanSummaryCDSc1sus2 DAC to ADC test

I ran a DAC to ADC test on c1sus2 channels where I hooked up the outputs on the DAC to the input channels on the ADC. We used different combinations of ADCs and DACs to make sure that there were no errors that cancel each other out in the end. I took a transfer function across these channel combinations to reproduce figure 1 in T2000188.

As seen in the two attached PDFs the channels seem to be working properly they have a flat response with a gain of 0.5 (-6 dB). This is the response that is expected and is the result of the DAC signal being sent as a single ended signal and the ADC receiving as a differential input signal. This should result in a recorded signal of 0.5 the amplitude of the actual output signal.

The drop off on the high frequency end is the result of the anti-aliasing filter and the anti-imaging filter. Both of these are 8-pole elliptical filters so when combined we should get a drop off of 320dB per decade. I measured the slope on the last few points of each filter and the averaged value was around 347dB per decade. This is slightly steeper than expected but since it is to cut off higher frequencies it shouldn't have an effect on the operation of the system. Also it is very close to the expected value.

The ripples seen before the drop off are also an effect of the elliptical filters and are seen in T2000188.

Note: the transfer function that doesn't seem to match the others is the heartbeat timing signal.

Attachment 1: data3_Plots.pdf
data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf
Attachment 2: data2_Plots.pdf
data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf
  16423   Fri Oct 22 17:35:08 2021 Ian MacMillanSummaryPEMParticle counter setup near BS Chamber

I have done some reading about where would be the best place to put the particle counter. The ISO standard (14644-1:2015) for cleanrooms is one every 1000 m^2 so one for every 30m x 30m space. We should have the particle counter reasonably close to the open chamber and all the manufactures that I read about suggest a little more than 1 every 30x30m. We will have it much closer than this so it is nice to know that it should still get a good reading. They also suggest keeping it in the open and not tucked away which is a little obvious. I think the best spot is attached to the cable tray that is right above the door to the control room. This should put it out of the way and within about 5m of where we are working. I ordered some cables to route it over there last night so when they come in I can put it up there. 

  16430   Tue Oct 26 18:24:00 2021 Ian MacMillanSummaryCDSc1sus2 DAC to ADC test

[Ian, Anchal, Paco]

After the Koji found that there was a problem with the power source Anchal and I fixed the power then reran the measurment. The only change this time around is that I increased the excitation amplitude to 100. In the first run the excitation amplitude was 1 which seemed to come out noise free but is too low to give a reliable value.

link to previous results

The new plots are attached.

Attachment 1: data2_Plots.pdf
data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf data2_Plots.pdf
Attachment 2: data3_Plots.pdf
data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf data3_Plots.pdf
  16447   Wed Nov 3 16:55:23 2021 Ian MacMillanSummarySUSSUS Plant Plan for New Optics

[Ian, Tega, Raj]

This is the rough plan for the testing of the new suspension models with the created plant model. We will test the suspensions on the plant model before we implement them into the full

  • Get State-space matrices from the surf model for the SOS. Set up simplant model on teststand
    • The state-space model is only 3 degrees of freedom. (even the surf's model)
    • There are filter modules that have the 6 degrees of freedom for the suspensions. We will use these instead. I have implemented them in the same suspension model that would hold the state space model. If we ever get the state space matrices then we can easily substitute them.
  • Load new controller model onto test stand. This new controller will be a copy of an existing suspension controller.
  • Hook up controller to simplant.  These should form a closed loop where the outputs from the controller go into the plant inputs and the plant outputs go to the controller inputs.
  • Do tests on set up.
    • Look at the step response for each degree of freedom. Then compare them to the results from an existing optic. 
      • Also, working with Raj let him do the same model in python then compare the two.
  • Make sure that the data is being written to the local frame handler.

MEDM file location

/opt/rtcds/userapps/release/sus/common/medm/hsss_tega_gautam

run using 

For ITMX display, use:

hsss_tega_gautam>medm -x -macro "site=caltech,ifo=c1,IFO=C1,OPTIC=ITMX,SUSTYPE=IM,DCU_ID=21,FEC=45" SUS_CUST_HSSS_OVERVIEW.adl

  16457   Mon Nov 8 17:52:22 2021 Ian MacMillanUpdateSUSSetting up suspension test model

[Ian, Tega]

We combined a controler and a plant model into a single modle (See first attachment) called x1sus_cp.mdl in the userapps folder of the cymac in c1sim. This model combines 2 blocks: the controler block which is used to control the current optics and is found in cvs/cds/rtcds/userapps/release/sus/c1/models/c1sus.mdl further the control block we are using comes from the same path but from the c1sup.mdl model. This plant model is the bases for all of my custom plant models and thus is a good starting point for the testing. It is also ideal because I know it can beeasily altered into a my state-space plant model. However, we had to make a few adjustments to get the model up to date for the cds system. So it is now a unique block.

These two library blocks are set in the userapps/lib folder on the cymac. This is the lib file that the docker system looks to when it is compiling models. For a quick overview see this. All other models have been removed from the MatLab path so that when we open x1sus_cp.mdl in MatLab it is using the same models it will compile with.

We could not find the rtbitget library part, but chris pointed us to userapps, and we copied it over using: scp /opt/rtcds/userapps/trunk/cds/common/models/rtbitget.mdl controls@c1sim:/home/controls/simLink/lib.

NOTE TO FUTURE IAN: don't forget that unit delays exist.

Next step: now that we have a model that is compiling and familiar we need to make medm screens. We will use the auto mdl2adl for this so that it is quick. Then we can start adding our custom pieces one by one so that we know that they are working. We will also work with Raj to get an independent python model working. Which will allow us to compare the cds and python models.

Attachment 1: x1sus_cp.png
x1sus_cp.png
  16460   Tue Nov 9 13:40:02 2021 Ian MacMillanSummaryComputer Scripts / ProgramsSUS Plant Plan for New Optics

[Ian, Tega]

After talking with Rana we have an updated plan. We will be working on this plan step by step in this order.

  1. Remove c1sim from the test stand rack and move it to the rack in the office next to the printer. When connecting it we will NOT connect it to the Martian network! This is to make sure that nothing is connected to the 40m system and we can't mess anything up.
  2. Once we have moved the computer over physically, we will need to update anyone who uses it on how to connect to it. The way we connect to it will have changed.
  3. Now that we have the computer moved and everyone can connect to it we will work on the model. Currently, we have the empty models connected.
    1. recompile the model since we moved the computer.
    2. verify that nothing has changed in the move and the model can still operate and compile properly
  4. The model has the proper structure but we need to fill it with the proper filters and such
    1. For the Plant model
      1. To get it up and running quickly we will use the premade plant filters for the plant model. These filters were made for the c1sup.mdl and should work in our modified plant model. This will allow us to verify that everything is working. And allow us to run tests on the system.
      2. We need to update the model and add the state space block. (we are skipping this step for now because we are fast-tracking the testing)  
        1. Check with Chris to make sure that this is the right way to do it. I am pretty sure it is, but I don't know anything
        2. Make the 6 DOF state-space matrix. We only have a three DOF one. The surf never made a 6 DOF. 
        3. Make the block to input into the model
        4. make a switch that will allow us to switch between the state-space model and the filter block
    2. For the controller
      1. Load filter coefficients for the controller model from one of the current optics and use this as a starting point.
      2. Add medm screens for the controller and plant. We are skipping this for now because we want results and we don't care if the screens look nice and are useable at the moment.
  5.  Test the model
    1. we will take an open-loop transfer function of all six of the DOFs to all other DOFs which will leave us with 36 TFs. Many will be zero
      1. If you are looking at this post then we are measuring transfer functions from the blue flags to the green flags across the plant model.
      2. We will want to look at the TFs across the controller

 

  16461   Tue Nov 9 16:55:52 2021 Ian MacMillanSummaryComputer Scripts / ProgramsSUS Plant Plan for New Optics

[Ian, Tega]

We have moved c1sim computer from the test stand to the server rack in the office area. (see picture)

It is connected to the general campus network. Through the network switch at the top of the rack. This switch seeds the entire Martian network.

Test to show that I am not lying:

  1. you can ping it or ssh into it at
    controls@131.215.114.116
    Using the same password as before. Notice this is not going through the nodus network.
  2. It also has a different beginning of the IP addresses. Martian network IP addresses start with 191.168.113

c1sim is now as connected to the 40m network as my mom's 10-year-old laptop.

unfortunately, I have not been able to get the x2go client to connect to it. I will have to investigate further. It is nice to have access to the GUI of c1sim occasionally.

Attachment 1: IMG_8107.JPG
IMG_8107.JPG
  16462   Tue Nov 9 18:05:03 2021 Ian MacMillanSummaryComputer Scripts / ProgramsSUS Plant Plan for New Optics

[Ian, Tega]

Now that the computer is in its new rack I have copied over the filter two files that I will use in the plant and the controller from pianosa:/opt/rtcds/caltech/c1/chans to the docker system in c1sim:/home/controls/docker-cymac/chans. That is to say, C1SUP.txt -> X1SUP.txt and C1SUS.txt -> X1SUS_CP.txt, where we have updated the names of the plant and controller inside the txt files to match our testing system, e.g. ITMX -> OPT_PLANT in plant model and ITMX -> OPT_CTRL in the controller and the remaining optics (BS, ITMY, PRM, SRM) are stripped out of C1SUS.txt in order to make X1SUS_CP.txt. 

Once the filter files were copied over need to add them to the filters that are in my models to do this I run the commands:

$  cd docker-cymac
$  eval $(./env_cymac)
$  ./login_cymac
 #  cd /opt/rtcds/tst/x1/medm/x1sus_cp
 #  medm -x X1SUS_OPT_PLANT_TM_RESP.adl

see this post for more detail

Unfortunately, the graphics forwarding from the docker is not working and is giving the errors:

arg: X1SUS_OPT_PLANT_TM_RESP.adl

locateResource 'X1SUS_OPT_PLANT_TM_RESP.adl'

isNetworkRequest X1SUS_OPT_PLANT_TM_RESP.adl

canAccess('X1SUS_OPT_PLANT_TM_RESP.adl', 4) = 0

can directly access 'X1SUS_OPT_PLANT_TM_RESP.adl'

isNetworkRequest X1SUS_OPT_PLANT_TM_RESP.adl

locateResource(X1SUS_OPT_PLANT_TM_RESP.adl...) returning 1

Error: Can't open display:

This means that the easiest way to add the filters to the model is through the GUI that can be opened through X2go client. It is probably easiest to get that working. graphics forwarding from inside the docker is most likely very hard. 

unfortunately again x2go client won't connect even with updated IP and routing. It gives me the error: unable to execute: startkde. Going into the files on c1sim:/usr/bin and trying to start startkde by myself also did not work, telling me that there was no such thing even though it was right in front of me.

  16466   Mon Nov 15 15:12:28 2021 Ian MacMillanSummaryComputer Scripts / ProgramsSUS Plant Plan for New Optics

[Ian, Tega]

We are working on three fronts for the suspension plant model:

  1. Filters
    1. We now have the state-space matrices as given at the end of this post. From these matrices, we can derive transfer functions that can be used as filter inputs. For a procedure see HERE. We accomplish this using Matlab's built-in ss(A,B,C,D); function. then we make it discrete using c2d(sys, 1/f); this gives us our discrete system running at the right frequency. We can get the transfer functions of either of these systems using tf(sys);
    2. from there we can copy the transfer functions into our photon filters.  Tega is working on this right now.
  2. State-Space
    1. We have our matrices as listed at the end of this post. With those compiled into a discrete system in MatLab we can use the code Chris made called rtss.m to convert this system into a .c file and a .h file.
    2. from there we have moved those files under the userapps folder in the docker system. then we added a c-code block to our .mdl model for the plant and pointed it at the custom c file we made. See section 7.2 of T080135-v10
    3. We have done all this and this should implement a custom state-space function into our .mdl file. the downside of this is that to change our SS model we have to edit the matrices we can't edit this from an medm screen. We have to recompile every time.
  3. Python Check
    1. This python check is run by Raj and will take in the state-space matrices which are given then will take transfer functions along all inputs and outputs and will compare them to what we have from the CDS model.

 

Here are the State-space matrices:

A=\begin{bmatrix} 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ -\omega_x^2(1+i/Q_{x}) & -\gamma_x & \omega_xb & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 \\ \frac{\omega_{\theta}^2}{l+b} & 0 & -\omega_{\theta}^2(1+i/Q_{\theta}) & -\gamma_{\theta} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & -\omega_{\phi}^2(1+i/Q_{\phi}) & -\gamma_{\phi} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \\ 0 & 0 & 0 & 0 & 0 & 0 & -\omega_{y}^2(1+i/Q_{y}) & -\gamma_{y}\end{bmatrix} 

  B=\begin{bmatrix} 0 & 0 & 0 & 0 \\ \frac{1}{m} & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & \frac{R_m}{I_{\theta}} & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & \frac{R_m}{I_{\phi}} & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{m} \end{bmatrix}      C=\begin{bmatrix} 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \end{bmatrix}      D=\begin{bmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{bmatrix}

A few notes: If you want the values for these parameters see the .yml file or the State-space model file. I also haven't been able to find what exactly this s is in the matrices.

UPDATE [11/16/21 4:26pm]: I updated the matrices to make them more general and eliminate the "s" that I couldn't identify. 

The input vector will take the form:

\begin{bmatrix} x \\ \dot{x} \\ \theta \\ \dot{\theta} \\ \phi \\ \dot{\phi} \\ y \\ \dot{y} \end{bmatrix}

where x is the position, theta is the pitch, phi is the yaw, and y is the y-direction displacement

 

 

  16469   Tue Nov 16 17:29:49 2021 Ian MacMillanSummaryComputer Scripts / ProgramsSUS Plant Plan for New Optics

[Ian, Tega]

Updated A, B, C, D matrices for the state-space model to remove bugs in the previous estimate of the system dynamics. Updated the last post to represent the current matrixes.

We used MatLab to get the correct time-series filter coefficients in ZPK format and added them to the filters running in the TM_RESP filter matrix.

Get the pos-pos transfer function from the CDS model. Strangely, this seems to take a lot longer than anticipated to generate the transfer function, even though we are mainly probing the low-frequency behavior of the system.  

For example, a test that should be taking approximately 6 minutes is taking well over an hour to complete. This swept sine (results below) was on the low settings to get a fast answer and it looks bad. This is a VERY basic system it shouldn't be taking this long to complete a Swept sine TF.

 

Noticed that we need to run eval $(./env_cymac) every time we open a new terminal otherwise CDS doesn't work as expected. Since this has been the source of quite a few errors already, we have decided to put it in the startup .bashrc script.  

loc=$(pwd)
cd ${HOME}/docker-cymac/
eval $(./env_cymac)
cd ${loc}
Attachment 1: x_x_TF1.pdf
x_x_TF1.pdf
  16477   Thu Nov 18 20:00:43 2021 Ian MacMillanSummaryComputer Scripts / ProgramsSUS Plant Plan for New Optics

[Ian, Raj, Tega]

Here is the comparison between the results of Raj's python model and the transfer function measurement done on the plant model by Tega and me.

As You can see in the graphs there are a few small spots of disagreement but it doesn't look too serious. Next we will measure the signals flowing through the entire plant and controller.

For a nicer (and printable) version of these plots look in the zipped folder under Plots/Plant_TF_Individuals.pdf

Attachment 1: Final_Plant_Testing.zip
  16481   Wed Nov 24 11:02:23 2021 Ian MacMillanSummaryComputersQuantization Noise Calculation Summary

I added mpmath to the quantization noise code. mpmath allows me to specify the precision that I am using in calculations. I added this to both the IIR filters and the State-space models although I am only looking at the IIR filters here. I hope to look at the state-space model soon. 

Notebook Summary:

I also added a new notebook which you can find HERE. This notebook creates a signal by summing two sine waves and windowing them.

Then that signal is passed through our filter that has been limited to a specific precision. In our case, we pass the same signal through a number of filters at different precisions.

Next, we take the output from the filter with the highest precision, because this one should have the lowest quantization noise by a significant margin, and we subtract the outputs of the lower precision filters from it. In summary, we are subtracting a clean signal from a noisy signal; because the underlying signal is the same, when we subtract them the only thing that should be left is noise. and since this system is purely digital and theoretical the limiting noise should be quantization noise.

Now we have a time series of the noise for each precision level (except for our highest precision level but that is because we are defining it as noiseless). From here we take a power spectrum of the result and plot it.

After this, we can calculate a frequency-dependent SNR and plot it. I also calculated values for the SNR at the frequencies of our two inputs. 

This is the procedure taken in the notebook and the results are shown below.

Analysis of Results:

The first thing we can see is that the precision levels 256 and 128 bits are not shown on our graph. the 256-bit signal was our clean signal so it was defined to have no noise so it cant be plotted. The 128-bit signal should have some quantization noise but I checked the output array and it contained all zeros. after further investigation, I found that the quantization noise was so small that when the result was being handed over from mpmath to the general python code it was rounding those numbers to zero. To overcome this issue I would have to keep the array as a mpmath object the entire time. I don't think this is useful because matplotlib probably couldn't handle it and it would be easier to just rewrite the code in C. 

The next thing to notice is sort of a sanity check thing. In general, low precision filters yield higher noise than high precision. This is a good quick sanity check. However, this does not hold true at the low end. we can see that 16-bit actually has the highest noise for most of the range. Chris pointed out that at low precisions that quantization noise can become so large that it is no longer a linearly coupled noise source. He also noted that this is prone to happen for low precision coefficients with features far below the Nyquist frequency like I have here. This is one explanation that seems to explain the data especially because this ambiguity is happening at 16-bit and lower as he points out. 

Another thing that I must mention, even if it is just a note to future readers, is that quantization noise is input dependent. by changing the input signal I see different degrees of quantization noise.

Analysis of SNR:

One of the things we hoped to accomplish in the original plan was to play around with the input and see how the results changed. I mainly looked at how the amplitude of the input signal scaled the SNR of the output. Below I include a table of the results. These results were taken from the SNR calculated at the first peak (see the last code block in the notebook) with the amplitude of the given sine wave given at the top of each column. this amplitude was given to both of the two sine waves even though only the first one was reported. To see an example, currently, the notebook is set up for measurement of input amplitude 10.

  0.1 Amplitude of input 1 Amplitude 100 Amplitude 1000 Amplitude
4-bit SNR 5.06e5 5.07e5 5.07e5 5.07e5
8-bit SNR 5.08e5 5.08e5 5.08e5 5.08e5
16-bit SNR 2.57e6 8.39e6 3.94e6 1.27e6
32-bit SNR 7.20e17 6.31e17 1.311e18 1.86e18
64-bit SNR 6.0e32 1.28e32 1.06e32 2.42e32
128-bit SNR unknown unknown unknown unknown

As we can see from the table above the SNR does not seem to relate to the amplitude of the input. in multiple instances, the SNR dips or peaks in the middle of our amplitude range.

 

Attachment 1: PSD_IIR_all.pdf
PSD_IIR_all.pdf
  16492   Tue Dec 7 10:55:25 2021 Ian MacMillanSummaryComputersQuantization Noise Calculation Summary

[Ian, Tega]

Tega and I have gone through the IIR Filter code and optimized it to make sure there aren't any areas that force high precision to be down-converted to low precision.

For the new biquad filter we have run into the issue where the gain of the filter is much higher than it should be. Looking at attachments 1 and 2, which are time series comparisons of the inputs and outputs from the different filters, we see that the scale for the output of the Direct form II filter shown in attachment 1 on the right is on the order of 10^-5 where the magnitude of the response of the biquad filter is on the order of 10^2. other than this gain the responses look to be the same. 

I am not entirely sure how this gain came into the system because we copied the c code that actually runs on the CDS system into python. There is a gain that affects the input of the biquad filter as shown on this slide of Matt Evans Slides. This gain, shown below as g, could decrease the input signal and thus fix the gain. However, I have not found any way to calculate this g.

 

 

With this gain problem we are left with the quantization noise shown in Attachment 4.

State Space:

I have controlled the state space filter to act with a given precision level. However, my code is not optimized. It works by putting the input state through the first state-space equation then integrating the result, which finally gets fed through the second state-space equation. 

This is not optimized and gives us the resulting quantization noise shown in attachment 5.

However, the state-space filter also has a gain problem where it is about 85 times the amplitude of the DF2 filter. Also since the state space is not operating in the most efficient way possible I decided to port the code chris made to run the state-space model to python. This code has a problem where it seems to be unstable. I will see if I can fix it

 

 

Attachment 1: DF2_TS.pdf
DF2_TS.pdf
Attachment 2: BIQ_TS.pdf
BIQ_TS.pdf
Attachment 4: PSD_COMP_BIQ_DF2.pdf
PSD_COMP_BIQ_DF2.pdf
Attachment 5: PSD_COMP_SS_DF2.pdf
PSD_COMP_SS_DF2.pdf
  16498   Fri Dec 10 13:02:47 2021 Ian MacMillanSummaryComputersQuantization Noise Calculation Summary

I am trying to replicate the simulation done by Matt Evans in his presentation  (see Attachment 1 for the slide in particular). 

He defines his input as x_{\mathrm{in}}=sin(2\pi t)+10^{-9} sin(2\pi t f_s/4) so he has two inputs one of amplitude 1 at 1 Hz and one of amplitude 10^-9 at 1/4th the sampling frequency  in this case: 4096 Hz

For his filter, he uses a fourth-order notch filter. To achieve this filter I cascaded two second-order notch filters (signal.iirnotch) both with locations at 1 Hz and quality factors of 1 and 1e6. as specified in slide 13 of his presentation

I used the same procedure outlined here. My results are posted below in attachment 2.

Analysis of results:

As we can see from the results posted below the results don't match. there are a few problems that I noticed that may give us some idea of what went wrong.

First, there is a peak in the noise around 35 Hz. this peak is not shown at all in Matt's results and may indicate that something is inconsistent.

the second thing is that there is no peak at 4096 Hz. This is clearly shown in Matt's slides and it is shown in the input spectrum so it is strange that it does not appear in the output.

My first thought was that the 4kHz signal was being entered at about 35Hz but even when you remove the 4kHz signal from the input it is still there. The spectrum of the input shown in Attachment 3 shows no features at ~35Hz.

The Input filter, Shown in attachment 4 shows the input filter, which also has no features at ~35Hz. Seeing how the input has no features at ~35Hz and the filter has no features at ~35Hz there must be either some sort of quantization noise feature there or more likely there is some sort of sampling effect or some effect of the calculation.

To figure out what is causing this I will continue to change things in the model until I find what is controlling it. 

I have included a Zip file that includes all the necessary files to recreate these plots and results.

Attachment 1: G0900928-v1_(dragged).pdf
G0900928-v1_(dragged).pdf
Attachment 2: PSD_COMP_BIQ_DF2.pdf
PSD_COMP_BIQ_DF2.pdf
Attachment 3: Input_PSD.pdf
Input_PSD.pdf
Attachment 4: Input_Filter.pdf
Input_Filter.pdf
Attachment 5: QuantizationN.zip
  2708   Wed Mar 24 12:38:17 2010 HartmutConfigurationGreen LockingBroadband PD for green PLL

Modified one of the PD assemblies carrying a large SI-Diode (~10mm diameter).

Removed elements used for resonant operation and changed PD readout to transimpedance

configuration. The opamp is a CLC409 with 240 Ohm feedback (i.e. transimpedance) resistor.

To prevent noise peaking at very high frequencies and get some decoupling of the PD,

I added a small series resistor in line with the PD and the inverting opamp input.

It was chosen as 13 Ohm, and still allows for operation up to ~100MHz.

Perhaps it could be smaller, but much more bandwith seems not possible with this opamp anyway.

Changes are marked in the schematic, and I list affected components here.

(Numbers refer to version 'PD327.SCH' from 30-April-1997):

-removed L4

-connected L3 (now open pad) via 100 Ohm to RF opamp output. This restores the DC sognal output.

-removed c17

-connected pin 3 of opamp via 25 Ohm to GND

-connected kathode of PD via 13 Ohm to pin 2 of opamp

-removed L6, C26, L5, C18, and C27

-shorted C27 pad to get signal to the RF output

 

Measured the optical TF with the test laser setup.

(Note that this is at 1064nm, although the PD is meant to work with green light at 532nm!)

Essentially it looks usable out to 100MHz, where the gain dropped only by about

6dB compared to 10MHz.

Beyond 100MHz the TF falls pretty steeply then, probably dominated by the opamp.

 

The maximal bias used is -150V.

If the bias is 'reduced' from -150V to -50V, the response goes down by 4dB at 10MHz and

by 9dB at 100MHz.

 The average output was 30mV at the RF output, corresponding to 60mV at the opamp output (50Ohm divider chain).

With 240 Ohm transimpedance this yields 250µA photo-current used for these transfer functions.

SiAmpl.png

 

SiPhase.png

 

 

  2745   Wed Mar 31 19:29:58 2010 HartmutUpdateElectronics(1cm-) Si PD transfer functions update

Recorded transfer functions for the 1cm Si-PD as described on p. 2708

for different biases. I put the plots in there, to keep the info in one place,

where the label on the PD case (which Steve made without asking him) points

to.

I talked to some people recently about the fact that the responsivity (A/W) of the PD

changes even at DC for different biases. I tested this again and should be more precise about this:

The first time I observed this was in the transfer functions as shown on p. 2708.

With 'DC' I meant 'low frequency' there, as you can still see an effect of the bias as low as 100kHz.

Then at one point I saw the responsivity changing with bias also at true DC.

However, it turned out that this is only the case if the photocurrent is too high.

If the photocurrent is 4mA, you need 400mV bias to get the max. responsivity.

For 2mA photocurrent, the responsivity is already maximal for 0V bias.

An effect for relative low frequencies remains however.

The DC check of responsivity was done with white light from a bulb.

 

 

  2752   Thu Apr 1 16:34:29 2010 HartmutUpdateGreen LockingSilicon PDs

just a few infos on Silicon PDs I looked up.

If you want to go beyond the 100MHz achievable with the device I worked on,

the one thing to improve is the opamp, where Steve is trying to find OPA657.

This is a FET with 1.6GHz BWP, minimum stable gain of 7, and 4.8nV/rt(Hz) noise.

Should be ok with 750-1000 Ohm transimpedance.

The other thing you might want to change is the PD

(although it might be the 1cm PD with high bias is as fast as smaller ones with lower bias).

There are two types of other Si diodes at the 40m right now (~3mm):

-Rana and I found a Centronic OSD 15-5T in the old equipment

-Frank gave me a Hamamatsu S1223-01 on a Thorlabs pre-amp device (could be taken out).

 

The Centronic OSD 15-5T has up to 80pF with 12 V bias according to the datasheet.

The Hamamatsu S1223-01 is stated with 20pF only, but stated to have a max. frequency resp. of 20MHz ('-3db point').

I dont know what this means, as the corner freq. of 10pF into 50Ohm is still 160MHz.

In any case there are faster 3mm types to start with, as for example Hamamatsu S3399 (~ 90$),

which is stated to have the corner at 100MHz with 50 Ohm load.

For this type the stated capacity (20pF) looks consistent with ~100MHz corner into 50 Ohm.

So probably you can get higher BW with this one using much smaller load, as in transimpedance stage.

 

 

  2757   Thu Apr 1 20:29:02 2010 HartmutUpdateGreen Lockingsimple PD test circuit

I made a simple PD test circuit which may allow to test PD response up to few 100MHz.

Its not for low noise, only for characterising PD response.

Here is the circuit:

The 2 capacitor values (for bypassing) are kind of arbitrary, just what I found around

(one medium, one small capacity). Could be improved by better RF types (e.g. Mica).

The PD type has no meaning. I put in the Centronic 15-T5 for a start.

The bias can be up to 20V for this diode.

The signal appears across R1. It is small, to make a large bandwidth.

R2 is just for slightly decoupling the signal from the following RF amplifier.

The wire into the RF amplifier is short (~cm). And the amplifier is supposed to have 50 Ohm

input impedance.

I use a mini circuits ZFL 500 here.

power supply for this is 15V.

pdtest.png

  10092   Tue Jun 24 13:04:49 2014 Harry, ManasaUpdateGeneralRazorblade Beam Analysis Setup

Harry will update this elog with details about his beam waist measurements for the old NPRO on the SP table.

  10099   Wed Jun 25 09:17:33 2014 Harry, ManasaUpdateGeneralRazorblade Beam Analysis Setup

Quote:

Harry will update this elog with details about his beam waist measurements for the old NPRO on the SP table.

 see http://nodus.ligo.caltech.edu:8080/40m/10098 for the update

  10083   Fri Jun 20 18:33:53 2014 HarryUpdateGeneralRazorblade Beam Analysis Setup

 Eric Q and I set up the optical configuration for razorblade beam analysis on SP table for future use.

It has been aligned, and will be in use on Monday.

The beam will be characterized for future characterization of optical fibers.

  10098   Wed Jun 25 09:16:52 2014 HarryUpdateGeneralRazorblade Measurements

Purpose

To use a razorblade to measure beam waist at multiple points along the optical axis, so as to later extrapolate the modal profile of the entire beam. This information will then be used to effectively couple AUX laser light to fibers for use in the frequency offset locking apparatus.

Data Acquisition

1) Step the micrometer-controlled razorblade across the beam at a given value of Z, along optical axis, in the plane orthogonal to it (arbitrarily called X).

2) At each value of X, record the corresponding output of a photodiode, (Thorlabs PD A55) here given in mV.

3) Repeat process at multiple points along Z

Analysis

Data from each iteration in the X were fitted to the error function shown below.

V(x) = A*(erf((x-m)/s)+c)

In the Y, they were fitted to:

V(x) = -A*(erf((x-m)/s)+c)

'A' corresponds to an amplitude, 'm' to a mean, 's' to a σ, and 'c' to an offset.

(Only because in Y measurements, the blade progressed toward eclipsing the beam, as opposed to in the X where it progressively revealed the beam.

These fits can be solved for x = (erf-1((V/A)-c)*s)+m1  which can be calculated at the points (Vmax/e2) and (Vmax*(1-1/e2)). The difference between these points will yield beam waist, w(z).

Conclusion

Calculations yielded waists of: X1=66.43um, X2=67.73um, X3=49.45um, Y1=61.20um, Y2=58.70, Y3=58.89

These data seem suspect, and shall be subjected to further analysis.

 

Attachment 1: 40m.zip
  10100   Wed Jun 25 09:30:44 2014 HarryUpdateGeneralWeekly Update

See attached weekly update

Attachment 1: Weekly_Update—June_25_thru_July_1.pdf
Weekly_Update—June_25_thru_July_1.pdf
  10103   Wed Jun 25 17:49:36 2014 HarryUpdateGeneralRazorblade Analysis Pt. 2

Reconfigured razorblade analysis setup on the PD table as per instructions. Used it to collect data to calculate beam waist with, analyses to follow.

See attached schematic for optical setup.

Attachment 1: RazorbladeSetup.pdf
RazorbladeSetup.pdf
  10106   Fri Jun 27 10:09:10 2014 HarryUpdateGeneralBeam Waist Measurement

Purpose

To use a razorblade to measure beam waist at four points along the optical axis, so as to later extrapolate the waist. This information will then be used to effectively couple AUX laser light to fibers for use in the frequency offset locking apparatus.

Data Acquisition

1) Step the micrometer-controlled razorblade across the beam at a given value of Z, along optical axis, in the plane orthogonal to it (arbitrarily called X).

2) At each value of X, record the corresponding output of a photodiode, (Thorlabs PD A55) here given in mV.

3) Repeat in Y plane at the same value of Z

4) Repeat process at multiple points along Z

Analysis

Data from each iteration were fitted to the error function shown below.

y(x) = (.5*P)*(1-erf((sqrt(2)*(x-x0))/wz))

'P' corresponds to peak power, 'x0' to the corresponding value of x (or y, as the case may be), and 'wz' to the spot size at the Z value in question.

The spot sizes from the four Z values were then fit to:

y(x) = w0*sqrt(1+((x*x)/(zr*zr)))

Where 'w0' corresponds to beam waist, and 'zr' to Rayleigh Range.

Conclusion

This yielded a Y-Waist of 783.5 um, and an X-Waist of 915.2 um.

The respective Rayleigh ranges were 2.965e+05 um (Y) and 3.145e+05  um (X).

Next

I will do the same analysis with light from the optical cables, which information I will then use to design a telescope to effectively couple the beams.

Attachment 1: BeamWaist.zip
  10139   Mon Jul 7 14:42:33 2014 HarryUpdateGeneralBeam Waists

I was finally able to get a reasonable measurement for the beam waist(s) of the spare NPRO.

Methods

I used a razorblade setup, pictured below, to characterize the beam waist of the spare 1064nm NPRO after a lens (PLCX-25.4-38.6-UV-1064) in order to subsequently calculate the overall waist of the beam. The setup is pictured below:

 

RzorbladeSetupFinal.pdf

After many failed attempts, this was the apparatus we (Manasa, Eric Q, Koji, and I) arrived with. The first lens after the laser was installed to focus the laser, because it's true waist was at an inaccessible location. Using the lens as the origin for the Z axis, I was able to determine the waist of the beam after the lens, and then calculate the beam waist of the laser itself using the equation wf = (lambda*f)/(pi*wo) where wf is the waist after the lens, lambda the wavelength of the laser, f the focal legth of the lens (75.0 mm in this case) and wo the waist before the lens.

We put the razorblade, second lens (to focus the beam onto the photodiode (Thorlabs PDA255)), and the PD with two attenuating filters with optical density of 1.0 and 3.0, all on a stage, so that they could be moved as a unit, in order to avoid errors caused by fringing effects caused by the razorblade.

I took measurements at six different locations along the optical axis, in orthogonal cross sections (referred to as X and Y) in case the beam turned to be elliptical, instead of perfectly circular in cross section. These measurements were carried out in 1" increments, starting at 2" from the lens, as measured by the holes in the optical table.

Analysis

Once I had the data, each cross section was fit to V(x) = (.5*Vmax)*(1-erf((sqrt(2)*(x-x0))/wz))+c, which corresponds to the voltage supplied to the PD at a particular location in x (or y, as the case may be). Vmax is the maximum voltage supplied, x0 is an offset in x from zero, wz is the spot size at that location in z, and c is a DC offset (ie the voltage on the PD when the laser is fully eclipsed.) These fits may all be viewed in the attached .zip file.

The spot sizes, extracted as parameters of the previous fits, were then fit to the equation which describes the propagation of the spot radius, w(z) = wo*sqrt(1+((z-b)/zr)^2)+c, w(z) = w0*sqrt(1+((((z-b)*.000001064)^2)/((pi*w0^2)^2))) where wo corresponds to beam waist, b is an offset in the z. Examples of these fits can be viewed in the attached .zip file.

Finally, since the waists given by the fits were the waists after a lens, I used the equation wf = (lambda*f)/(pi*wo), described above, to determine the waist of the beam before the lens.

Plots

 

note: I was not able to open the first measurement in the X plane (Z = 2in). The rest of the plots have been included in the body of the elog, as per Manasa's request.

Y_Waist_Fit.pngY_plane_measurement_6_(Z_7in).pngY_plane_measurement_5_(Z_6in).pngY_plane_measurement_4(Z_5in).pngY_plane_measurement_3_(Z_4in).pngY_plane_measurement_2_(Z_3in).png

 

Y_plane_measurement_1_(Z_2in).png X_plane_measurement_6_(Z_7in).pngX_plane_measurement_5_(Z_6in).pngX_plane_measurement_4_(Z_5in).pngX_plane_measurement_3_(Z_4in).pngX_plane_measurement_2_(Z_3in).pngX_plane_measurement_1_(Z_2in).pngX_Waist_fit.png

Conclusion

The X Waist after the lens (originally yielded from fit parameters) was 90.8 27.99 ± .14 um. The corresponding Y Waist was 106.2 30.22 ± .11 um.

After adjustment for the lens, the X Waist was 279.7 907.5 ± 4.5 um and the Y Waist was 239.2 840.5 ± 3.0 um.

edit: After making changes suggested by koji, these were the new results of the fits.

Attachments

Attached you should be able to find the razor blade schematic, all of the fits, along with code used to generate them, plus the matlab workspace containing all the necessary variables.

NOTE: Rana brought to my attention that my error bars need to be adjusted, which I will do as soon as possible.

Attachment 2: erFitFinal2.zip
  10158   Tue Jul 8 23:59:49 2014 HarryUpdateGeneralWeekly Plan (7.8.14)

 Last Week:

-I continued to struggle with the razorblade beam analysis, though after a sixth round of measurements, and a lot of fiddling around with fit parameters in matlab, there seems to be a light at the end of the tunnel.

 

Next Week:

-I plan to check my work with the beamscan tomorrow (wednesday) morning

-Further characterize the light from the fibers, and set up the collimator

-Design and hopefully construct the telescope that will focus the beam into the collimator

 

Materials:

- Razorblade setup or beamscan (preferably beamscan)

- Fiber Illuminator

- Collimator (soon to be ordered)

- Lenses for telescope (TBD)

 

  10170   Wed Jul 9 23:02:33 2014 HarryUpdateGeneralBeam Waists via Beamscan

 Today, I borrowed the beam profiler from Brian (another SURF) in order to double check my razor blade measurement figures, using the below setup:

beamscanSetup.png

Measurements are included in the a la mode code that is attached entitled beamfit.m. The beam fitting application yielded me waists (after the lens) of 35.44 um in the x plane, and 33.26 um in the y plane. These are both within 3 um of the measurements I found using the razor blade method. (I moved and resized the labels for the waists in the figure below for readability purposes.)

beamScanWaists.png

I then plugged these waists back into ALM, in addition to the lens specifications, to determine waist size and location of the NPRO, which turned out to be 543 um in the x located at Z = 1.160m, and 536 um in the y, located at 1.268m. These measurements are based upon zero at the waist after the lens, and the positive direction being back toward the NPRO.

beamProfileX.pngbeamProfileY.png

The only systemic difference between these measurement and my original razor blade measurements was that I had taken the focal length of the lens as 75mm, which is advertised on the manufacturer's site. However, the more detailed specs revealed that the focal length was 85.8mm at 1064nm, which made a difference of about 400 um for the final waist determination.

Attachment 4: beamScanWaist.zip
  10175   Thu Jul 10 15:27:09 2014 HarryUpdateGeneralCoupling telescope design

 I designed this telescope to couple the 1064 NPRO into the PM980 fiber, using lenses from the Thorlabs LSB04-C kit.

The collimator is a CFC-2X-C, which has a variable focus length (2.0, 4.6, 7.5, and 11.0 mm) which gives corresponding angles of divergence of 0.298, 0.130, 0.79, and 0.054 degrees by the formula theta = (180*MFD) / (pi*f).

Then, using these values I calculated the spot size of a beam collimated by the CFC-2X-C, using f = w / tan(theta) where w is the spot size. This gave a value of 10.4 um.

I used this value (10.4 um) as a target waist for the telescope system, with the NPRO waist as a seed, at the origin.

telescopeBeamProfile.png

It consists of two lenses, one located at Z = 77cm f = 50cm, and the second located at Z = 85.88 cm f = 2.54cm, which yields a waist of 13um at Z = 88.32cm, (which is where the collimator would go) for an overlap of .974.

Note that the telescope is so far "downrange" from the NPRO waist because it's a virtual waist, and the NPRO itself is located at about Z = 73cm.

Find attached the alm code used.

Attachment 2: telescope.zip
  10185   Fri Jul 11 15:29:26 2014 HarryUpdateGeneralTelescope With Collimator

 I used a la mode to make a design for the coupling telescope with a 3.3um target waist, that included the collimator in the overall design. The plot is below, and code is attached.

The components are as follows:

 label         z (m)     type    parameters         

   -----         -----     ----    ----------         

    lens1         0.7681    lens    focalLength: 0.5000

    lens2         0.8588    lens    focalLength: 0.0350

    collimator    0.8832    lens    focalLength: 0.0020

 

The z coordinates are as measured from the beam waist of the NPRO (the figure on the far left  of the plot).

Moving forward, this setup will be used to couple the NPRO (more specifically, the AUX lasers) light into the SM 980 fibers, as well as to help characterize the fibers themselves.

Ultimately, this will be a key component in the Frequency Offset Locking project that Akhil and I are working on, as it will transport the AUX light to the PSL, where the two beams will be beaten with each other to generate the input signal to the PID control loop, which will actuate the temperature servos of the AUX lasers.

Attachment 1: telescope_2.zip
Attachment 2: telescopeWCollimator.png
telescopeWCollimator.png
ELOG V3.1.3-