Merge pull request #12 from jacotay7/dev

Huge merge of documentation and features from the past few months.

.github/workflows/draft-pdf.yml

@@ -13,7 +13,7 @@ jobs:
           journal: joss
           paper-path: paper.md
       - name: Upload
-        uses: actions/upload-artifact@v1
+        uses: actions/upload-artifact@v4
         with:
           name: paper
           path: paper.pdf

README.md

@@ -34,8 +34,49 @@ cd pyRTC
 pip install .
 ```
 
+Optionally, install with docs libraries:
+
+```
+cd pyRTC
+pip install .[docs]
+```
+
 # Getting Started 
 
+
+## Documentation
+
+You can build the docs locally by following these steps:
+
+First, install the required packages, this can be done by:
+
+```
+cd pyRTC
+pip install .[docs]
+```
+
+See Installation Instructions for more info.
+
+Next, navigate to the pyRTC folder and build the docs:
+
+```
+cd ~/pyRTC/docs/source
+make html
+```
+
+The docs will be built in the `docs/source/_build/html` folder. You can run them locally with:
+
+```
+sphinx-autobuild . _build/html
+```
+
+This process can also be done using the including `build_and_run.sh` script located in the `docs/source` folder.
+
+```
+cd ~/pyRTC/docs/source
+./build_and_run.sh
+```
+
 ## Hardware
 
 Since each AO system utilizes a unique hardware configuration, every new AO system will require a programmer to write interfacing code to control the hardware. This is unavoidable. Therefore, one of pyRTC's goals is compile a fairly complete set of API implementations for common hardware components like ALPAO DM's, PI modulators, and FLIR cameras so that community members will have easy access to existing exemplars. However, there is a significantly prohibitive cost associated with buying one of each possible hardware component an AO research might need to interact with. Therefore, we will rely on the community to provide their pyRTC implementations for their hardware components so that the whole community can benefit.
@@ -52,7 +93,6 @@ For our examples we will be using the open-source AO simulation software OOPAO d
 
 We have provided an example of how to set up a single conjugate AO (SCAO) system using a Pyramid Wavefront Sensor. This example uses the OOPAO simulation software to simulate the hardware components of an AO system so that pyRTC can be tested without the need of real AO hardware. All examples can be found under the `examples` folder.
 
-
 # Contributing
 
 pyRTC is an open-source software intended to be built for and by the AO community. Please use pyRTC freely and create new branches for your systems. We ask that if you implement new features that may be of broader interest to the community while integrating pyRTC into your AO system, please contribute to the code base by opening up an issue or pull request on github.   

SHARP_LAB/RTC.py

@@ -2,15 +2,17 @@
 #Import pyRTC classes
 from pyRTC.Pipeline import *
 from pyRTC.utils import *
+from pyRTC.hardware import *
+import matplotlib.pyplot as plt
 import os
 os.chdir("/home/whetstone/pyRTC/SHARP_LAB")
 RECALIBRATE = False
-
+CLEAR_SHMS =  False
 # %% Clear SHMs
-# from pyRTC.Pipeline import clear_shms
-# shm_names = ["signal"]
-# shm_names = ["wfs", "wfsRaw", "wfc", "wfc2D", "signal", "signal2D", "psfShort", "psfLong"] #list of SHMs to reset
-# clear_shms(shm_names)
+if CLEAR_SHMS:
+    from pyRTC.Pipeline import clear_shms
+    shm_names = ["wfs", "wfsRaw", "wfc", "wfc2D", "signal", "signal2D", "psfShort", "psfLong"] #list of SHMs to reset
+    clear_shms(shm_names)
 
 # %% IMPORTS
 # config = '/home/whetstone/pyRTC/SHARP_LAB/config_SR.yaml'
@@ -38,14 +40,13 @@
 loop.launch()
 
 # %% OPTIMIZERS
-from pyRTC.hardware.PIDOptimizer import PIDOptimizer
 pidOptim = PIDOptimizer(read_yaml_file(config)["optimizer"]["pid"], loop)
-
-from pyRTC.hardware.NCPAOptimizer import NCPAOptimizer
+loopOptim = loopOptimizer(read_yaml_file(config)["optimizer"]["loop"], loop)
 ncpaOptim = NCPAOptimizer(read_yaml_file(config)["optimizer"]["ncpa"], loop, slopes)
 
 # %% Calibrate
 
+
 if RECALIBRATE == True:
 
     slopes.setProperty("refSlopesFile", "")
@@ -62,7 +63,7 @@
     input("Sources On?")
     input("Is Atmosphere Out?")
 
-    slopes.run("computeImageNoise")
+    # slopes.run("computeImageNoise")
     slopes.run("takeRefSlopes")
     slopes.setProperty("refSlopesFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/ref.npy")
     slopes.run("saveRefSlopes")
@@ -83,40 +84,43 @@
     wfc.run("flatten")
     time.sleep(1)
 
-    # input("Is Atmosphere In?")
-    # #  DOCRIME OL
-    # loop.setProperty("IMMethod", "docrime")
-    # loop.setProperty("delay", 1)
-    # loop.setProperty("pokeAmp", 2e-2)
-    # loop.setProperty("numItersIM", 30000)
-    # loop.setProperty("IMFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/docrime_IM.npy")
-    # wfc.run("flatten")
-    # loop.run("computeIM")
-    # loop.run("saveIM")
-    # wfc.run("flatten")
-    # time.sleep(1)
+    input("Is Atmosphere In?")
+    #  DOCRIME OL
+    loop.setProperty("IMMethod", "docrime")
+    loop.setProperty("delay", 2) #Needs to be set in the CONFIG
+    loop.setProperty("pokeAmp", 1e-2)
+    loop.setProperty("numItersIM", 1000)
+    loop.setProperty("IMFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/IM_OL_docrime.npy")
+    wfc.run("flatten")
+    loop.run("computeIM")
+    loop.run("saveIM")
+    wfc.run("flatten")
+    time.sleep(1)
 
 
 
 # %% Compute CM
 loop.setProperty("IMFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/IM.npy")
+# loop.setProperty("IMFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/OL_DOCRIME.npy")
+# loop.setProperty("IMFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/OL_DOCRIME_CL_docrime.npy")
+# loop.setProperty("IMFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/ESCAPE.npy")
+loop.setProperty("numDroppedModes", 0)
+loop.setProperty("gain",0.40)
+loop.setProperty("leakyGain", 0.021)
 loop.run("loadIM")
 time.sleep(0.5)
-loop.setProperty("numDroppedModes", 15)
-loop.run("computeCM")
-time.sleep(0.5)
 # %% Adjust Gains
-# loop.run("setGain",1e-2)
-# loop.setProperty("leakyGain", 0.3)
-loop.setProperty("leakyGain", 1e-2)
+# loop.run("setGain",3e-1)
+# # loop.setProperty("leakyGain", 0.3)
+# loop.setProperty("leakyGain", 1e-2)
 #PID GAINS (only for PID integrator)
-loop.setProperty("pGain", 0.3)
-loop.setProperty("iGain", 5e-2)
-loop.setProperty("dGain", 5e-2)
-loop.setProperty("controlLimits", [-0.05, 0.05])
-loop.setProperty("absoluteLimits", [-1.0, 1.0])
-loop.setProperty("integralLimits", [-1.0, 1.0])
-# %%Launch Loop for 5 seconds
+# loop.setProperty("pGain", 0.3)
+# loop.setProperty("iGain", 5e-2)
+# # loop.setProperty("dGain", 5e-2)
+# loop.setProperty("controlLimits", [-0.05, 0.05])
+# loop.setProperty("absoluteLimits", [-0.3, 0.3])
+# loop.setProperty("integralLimits", [-0.1, 0.1])
+# %%Launch Loop for 5 seconds ###########################################################
 wfc.run("flatten")
 loop.run("start")
 time.sleep(1)
@@ -125,8 +129,58 @@
 wfc.run("flatten")
 wfc.run("flatten")
 
+
+# %% CL DOCRIME Start
+loop.setProperty("pokeAmp", 0.003)
+loop.setProperty("clDocrime", True)
+
+# %% SOLVE CL DOCRIME
+loop.run("solveDocrime")
+
+# %% CL DOCRIME Stop
+loop.setProperty("clDocrime", False)
+
+# %% Plot DOCRIME
+
+
+baseline = np.load("../SHARP_LAB/calib/baseline.npy")
+OLDOCRIME = np.load("../SHARP_LAB/calib/OL_DOCRIME.npy")
+CLDOCRIME = np.load("../SHARP_LAB/calib/CL_DOCRIME.npy")
+ESCAPE = np.load("../SHARP_LAB/calib/ESCAPE.npy")
+im = ESCAPE
+plt.imshow(im, aspect="auto")
+plt.show()
+
+a,b,c = np.std(baseline,axis = (0)), np.std(OLDOCRIME,axis = (0)), np.std(CLDOCRIME,axis = (0))
+d = np.std(ESCAPE,axis = (0))
+
+plt.plot(a, label = 'baseline')
+plt.plot(b, label = 'OL DOCRIME')
+plt.plot(c, label = 'CL DOCRIME')
+plt.plot(d, label = 'ESCAPE')
+plt.legend()
+plt.show()
+
+plt.title("Modal Gains")
+plt.plot(b/a, label = 'OL DOCRIME')
+plt.plot(c/a, label = 'CL DOCRIME')
+plt.plot(d/a, label = 'ESCAPE')
+plt.legend()
+plt.show()
+
+vsa = np.load("../SHARP_LAB/calib/validSubAps.npy")
+
+imFull = np.zeros((im.shape[-1],*vsa.shape))
+for i in range(imFull.shape[0]):
+    imFull[i][vsa] = im[:,i]
+N = 0
+for i in range(N,N + 5):
+    plt.imshow(imFull[i])
+    plt.colorbar()
+    plt.show()
+
 #%% Optimize PID
-pidOptim.numReads = 10
+pidOptim.numReads = 20
 pidOptim.numSteps = 50
 pidOptim.isPOL = False
 pidOptim.mode = "tiptilt"
@@ -137,20 +191,49 @@
     pidOptim.optimize()
 pidOptim.applyOptimum()
 #%% Optimize NCPA
-ncpaOptim.numReads = 2
-ncpaOptim.startMode = 2
-ncpaOptim.endMode = 30
-ncpaOptim.optimize()
-# ncpaOptim.applyOptimum()
-# wfc.run("saveShape")
-# slopes.run("takeRefSlopes")
-# slopes.setProperty("refSlopesFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/ref.npy")
-# slopes.run("saveRefSlopes")
-# psfCam.run("takeModelPSF")
-# psfCam.setProperty("modelFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/modelPSF.npy")
-# psfCam.run("saveModelPSF")
-
+#%
+import optuna
+optuna.logging.set_verbosity(optuna.logging.WARNING)
+numOptim = 3
+maxAMP = 0.02
+amps = np.linspace(maxAMP, maxAMP/5, numOptim)
+for i in range(numOptim):
+    ncpaOptim.resetStudy()
+    psfCam.setProperty("integrationLength", 5)
+    time.sleep(2)
+    ncpaOptim.numReads = 3
+    ncpaOptim.startMode = 0
+    ncpaOptim.endMode = 20 #wfc.getProperty("numModes")
+    ncpaOptim.numSteps = 1500
+    ncpaOptim.correctionMag = amps[i]
+    ncpaOptim.isCL = False
+    for i in range(1):
+        ncpaOptim.optimize()
+    ncpaOptim.applyNext()
+
+    wfc.run("saveShape")
+    # slopes.run("takeRefSlopes")
+    # slopes.setProperty("refSlopesFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/ref.npy")
+    # slopes.run("saveRefSlopes")
+    psfCam.setProperty("integrationLength", 2000)
+    time.sleep(2)
+    psfCam.run("takeModelPSF")
+    psfCam.setProperty("modelFile", "/home/whetstone/pyRTC/SHARP_LAB/calib/modelPSF_PyWFS.npy")
+    psfCam.run("saveModelPSF")
+    wfc.run("loadFlat")
+
 
+#%% Loop Optimizer
+psfCam.setProperty("integrationLength", 100)
+time.sleep(1)
+loopOptim.numReads = 10
+# loopOptim.maxGain = 0.6
+# loopOptim.maxLeak = 0.1
+# loopOptim.maxDroppedModes = 40
+loopOptim.numSteps = 300
+for i in range(1):
+    loopOptim.optimize()
+loopOptim.applyOptimum()
 
 
 # %% Plots
@@ -262,14 +345,14 @@
 # %%
 from tqdm import tqdm
 
-folder = "~/Downloads/robin-april-16/"
+folder = "~/Downloads/robin-aug/"
 
 # numModes = 10
 # startMode = 0
 # endMode = wfc.numModes - 1
 filelist = ['cnnx2_phase.npy', 'cnnx4_phase', 'linx2_phase.npy', 'linx4_phase.npy']
-N = 4
-numModes = 11
+N = 1
+numModes = 3
 RANGE = 2
 modelist = np.linspace(-RANGE, RANGE, numModes) #.astype(int)
 
@@ -280,6 +363,8 @@
     wfc.flatten()
     time.sleep(0.1)
     d = np.read(f'{folder}/{ff}')
+
+
     for i, mode in enumerate(modelist): #range(numModes):
         correction = np.zeros_like(wfc.read())
         for j in range(N):
@@ -293,8 +378,8 @@
             wfc.flatten()
             time.sleep(0.1)
 
-    np.save(f'{folder}/psfs_{ff}', psfs)
-    np.save(f'{folder}/cmds_{ff}', cmds)
+    np.savez(f'{folder}/tiptiltcalib_{ff}', psfs, cmds)
+    # np.save(f'{folder}/cmds_{ff}', cmds)
 # %%
 # %% CODE TO UPDATE FLAKY SUB APERTURES
 
@@ -337,4 +422,48 @@
 # %%
 for trial in pidOptim.study.trials:
     print(trial)
-    break
+    break
+
+# %% Save long PSF
+
+##CHANGE THIS LINE
+filename = "BASLINE_LONG_PSF"
+numReads = 100
+folder = "/home/whetstone/pyRTC/SHARP_LAB/data/"
+filename = folder + filename
+
+
+psfShm, _, _ = initExistingShm("psfLong")
+strehlShm, _, _ = initExistingShm("strehl")
+mean_psf = np.zeros_like(psfShm.read_noblock())
+mean_strehl = 0
+
+for i in range(numReads):
+    mean_psf += psfShm.read()
+    mean_strehl += strehlShm.read()[0]
+
+mean_psf /= numReads
+mean_strehl /= numReads
+mean_strehl = np.round(mean_strehl, 2)
+filename += str(mean_strehl)
+
+np.save(filename,mean_psf)
+
+plt.imshow(mean_psf)
+plt.title(f"SR: {mean_strehl}")
+plt.show()
+# %%
+
+# %%
+
+
+shm1, _, _ = initExistingShm("wfsRaw")
+N = 100000
+start = time.time()
+for i in range(N):
+    precise_delay(0.1)
+    shm1.read_noblock(SAFE=False)
+execTime = 1e6*(time.time()-start)/N
+
+print(f"Mean Execution Time: {execTime}us")
+# %%