|
| 1 | +{ |
| 2 | + "cells": [ |
| 3 | + { |
| 4 | + "cell_type": "markdown", |
| 5 | + "id": "dc42d337-8478-43d2-a3c9-562d6ccbe5b6", |
| 6 | + "metadata": {}, |
| 7 | + "source": [ |
| 8 | + "# CLT-006: Control Loop Single and Double Exposure Test\n", |
| 9 | + "\n", |
| 10 | + "Owner: **Bryce Kalmbach** <br>\n", |
| 11 | + "Last Verified to Run: **2025-04-03** <br>\n", |
| 12 | + "Software Version:\n", |
| 13 | + " - `ts_wep`: **14.1.1**\n", |
| 14 | + " - `donut_viz`: **1.6.2**\n", |
| 15 | + " - `lsst_distrib`: **w_2025_13**\n", |
| 16 | + "\n", |
| 17 | + "## Test Details:\n", |
| 18 | + "In this experiment we calculate the closed loop with 15 second exposure then reset the state back to the initial state and re-run the closed loop with 30 seconds.\n", |
| 19 | + "We then compare the residual AOS FWHM during the progression of the closed loop." |
| 20 | + ] |
| 21 | + }, |
| 22 | + { |
| 23 | + "cell_type": "code", |
| 24 | + "execution_count": null, |
| 25 | + "id": "a5e0fbcf", |
| 26 | + "metadata": {}, |
| 27 | + "outputs": [], |
| 28 | + "source": [ |
| 29 | + "# Times Square Parameters\n", |
| 30 | + "collection_name = 'u/brycek/aosRefitWcs_danish_singleBlends_80pxMinSep'\n", |
| 31 | + "day_obs = 20241118\n", |
| 32 | + "min_seq_num_15_sec = 31\n", |
| 33 | + "max_seq_num_15_sec = 53\n", |
| 34 | + "min_seq_num_30_sec = 54\n", |
| 35 | + "max_seq_num_30_sec = 75" |
| 36 | + ] |
| 37 | + }, |
| 38 | + { |
| 39 | + "cell_type": "code", |
| 40 | + "execution_count": null, |
| 41 | + "id": "8f7a6e94-7e2f-4a68-a56a-0b6a6fe1e3e7", |
| 42 | + "metadata": {}, |
| 43 | + "outputs": [], |
| 44 | + "source": [ |
| 45 | + "import numpy as np\n", |
| 46 | + "from copy import copy\n", |
| 47 | + "from matplotlib import pyplot as plt\n", |
| 48 | + "from lsst.daf.butler import Butler\n", |
| 49 | + "\n", |
| 50 | + "# Can uncomment when running outside Times Square (Hopefully temporary)\n", |
| 51 | + "#from lsst.ts.wep.utils import getPsfGradPerZernike\n", |
| 52 | + "\n", |
| 53 | + "from astropy.table import Table\n", |
| 54 | + "%matplotlib inline" |
| 55 | + ] |
| 56 | + }, |
| 57 | + { |
| 58 | + "cell_type": "code", |
| 59 | + "execution_count": null, |
| 60 | + "id": "babf357a-875a-4766-ad88-8a3e61e0147a", |
| 61 | + "metadata": {}, |
| 62 | + "outputs": [], |
| 63 | + "source": [ |
| 64 | + "butler = Butler('/repo/main')" |
| 65 | + ] |
| 66 | + }, |
| 67 | + { |
| 68 | + "cell_type": "code", |
| 69 | + "execution_count": null, |
| 70 | + "id": "9878cab0-402f-4ac0-a6c9-c7852ac3ff0b", |
| 71 | + "metadata": {}, |
| 72 | + "outputs": [], |
| 73 | + "source": [ |
| 74 | + "camera = butler.get('camera', {'instrument': \"LSSTComCam\"}, collections=collection_name)" |
| 75 | + ] |
| 76 | + }, |
| 77 | + { |
| 78 | + "cell_type": "markdown", |
| 79 | + "id": "c1b56ecb-ce5c-4162-89d4-abf040b309aa", |
| 80 | + "metadata": {}, |
| 81 | + "source": [ |
| 82 | + "## Define Functions Needed (Temporary)\n", |
| 83 | + "This is hopefully temporary since this function lives in `ts_wep` which is currently unavailable from Times Square" |
| 84 | + ] |
| 85 | + }, |
| 86 | + { |
| 87 | + "cell_type": "code", |
| 88 | + "execution_count": null, |
| 89 | + "id": "35405e6c-5b07-49d2-85c9-2dfd4f3bd089", |
| 90 | + "metadata": {}, |
| 91 | + "outputs": [], |
| 92 | + "source": [ |
| 93 | + "import galsim\n", |
| 94 | + "\n", |
| 95 | + "def getPsfGradPerZernike(\n", |
| 96 | + " diameter: float = 8.36,\n", |
| 97 | + " obscuration: float = 0.612,\n", |
| 98 | + " jmin: int = 4,\n", |
| 99 | + " jmax: int = 22,\n", |
| 100 | + ") -> np.ndarray:\n", |
| 101 | + " \"\"\"Get the gradient of the PSF FWHM with respect to each Zernike.\n", |
| 102 | + "\n", |
| 103 | + " This function takes no positional arguments. All parameters must be passed\n", |
| 104 | + " by name (see the list of parameters below).\n", |
| 105 | + "\n", |
| 106 | + " Parameters\n", |
| 107 | + " ----------\n", |
| 108 | + " diameter : float, optional\n", |
| 109 | + " The diameter of the telescope aperture, in meters.\n", |
| 110 | + " (the default, 8.36, corresponds to the LSST primary mirror)\n", |
| 111 | + " obscuration : float, optional\n", |
| 112 | + " Central obscuration of telescope aperture (i.e. R_outer / R_inner).\n", |
| 113 | + " (the default, 0.612, corresponds to the LSST primary mirror)\n", |
| 114 | + " jmin : int, optional\n", |
| 115 | + " The minimum Noll index, inclusive. Must be >= 0. (the default is 4)\n", |
| 116 | + " jmax : int, optional\n", |
| 117 | + " The max Zernike Noll index, inclusive. Must be >= jmin.\n", |
| 118 | + " (the default is 22.)\n", |
| 119 | + "\n", |
| 120 | + " Returns\n", |
| 121 | + " -------\n", |
| 122 | + " np.ndarray\n", |
| 123 | + " Gradient of the PSF FWHM with respect to the corresponding Zernike.\n", |
| 124 | + " Units are arcsec / micron.\n", |
| 125 | + "\n", |
| 126 | + " Raises\n", |
| 127 | + " ------\n", |
| 128 | + " ValueError\n", |
| 129 | + " If jmin is negative or jmax is less than jmin\n", |
| 130 | + " \"\"\"\n", |
| 131 | + " # Check jmin and jmax\n", |
| 132 | + " if jmin < 0:\n", |
| 133 | + " raise ValueError(\"jmin cannot be negative.\")\n", |
| 134 | + " if jmax < jmin:\n", |
| 135 | + " raise ValueError(\"jmax must be greater than jmin.\")\n", |
| 136 | + "\n", |
| 137 | + " # Calculate the conversion factors\n", |
| 138 | + " conversion_factors = np.zeros(jmax + 1)\n", |
| 139 | + " for i in range(jmin, jmax + 1):\n", |
| 140 | + " # Set coefficients for this Noll index: coefs = [0, 0, ..., 1]\n", |
| 141 | + " # Note the first coefficient is Noll index 0, which does not exist and\n", |
| 142 | + " # is therefore always ignored by galsim\n", |
| 143 | + " coefs = [0] * i + [1]\n", |
| 144 | + "\n", |
| 145 | + " # Create the Zernike polynomial with these coefficients\n", |
| 146 | + " R_outer = diameter / 2\n", |
| 147 | + " R_inner = R_outer * obscuration\n", |
| 148 | + " Z = galsim.zernike.Zernike(coefs, R_outer=R_outer, R_inner=R_inner)\n", |
| 149 | + "\n", |
| 150 | + " # We can calculate the size of the PSF from the RMS of the gradient of\n", |
| 151 | + " # the wavefront. The gradient of the wavefront perturbs photon paths.\n", |
| 152 | + " # The RMS quantifies the size of the collective perturbation.\n", |
| 153 | + " # If we expand the wavefront gradient in another series of Zernike\n", |
| 154 | + " # polynomials, we can exploit the orthonormality of the Zernikes to\n", |
| 155 | + " # calculate the RMS from the Zernike coefficients.\n", |
| 156 | + " rms_tilt = np.sqrt(np.sum(Z.gradX.coef**2 + Z.gradY.coef**2) / 2)\n", |
| 157 | + "\n", |
| 158 | + " # Convert to arcsec per micron\n", |
| 159 | + " rms_tilt = np.rad2deg(rms_tilt * 1e-6) * 3600\n", |
| 160 | + "\n", |
| 161 | + " # Convert rms -> fwhm\n", |
| 162 | + " fwhm_tilt = 2 * np.sqrt(2 * np.log(2)) * rms_tilt\n", |
| 163 | + "\n", |
| 164 | + " # Save this conversion factor\n", |
| 165 | + " conversion_factors[i] = fwhm_tilt\n", |
| 166 | + "\n", |
| 167 | + " return conversion_factors[jmin:]\n" |
| 168 | + ] |
| 169 | + }, |
| 170 | + { |
| 171 | + "cell_type": "markdown", |
| 172 | + "id": "15c10091-61d9-4e79-acbf-39aa16e68174", |
| 173 | + "metadata": {}, |
| 174 | + "source": [ |
| 175 | + "## Gather Visit Data" |
| 176 | + ] |
| 177 | + }, |
| 178 | + { |
| 179 | + "cell_type": "code", |
| 180 | + "execution_count": null, |
| 181 | + "id": "efed4f97-5e8a-4c93-93a0-d5124368a658", |
| 182 | + "metadata": {}, |
| 183 | + "outputs": [], |
| 184 | + "source": [ |
| 185 | + "visit_tables_15_sec = butler.query_datasets('aggregateAOSVisitTableAvg', \n", |
| 186 | + " collections=collection_name,\n", |
| 187 | + " where=f\"exposure.day_obs = {day_obs} and exposure.seq_num >= {min_seq_num_15_sec} and exposure.seq_num <= {max_seq_num_15_sec} and instrument = 'LSSTComCam'\")\n", |
| 188 | + "visit_tables_30_sec = butler.query_datasets('aggregateAOSVisitTableAvg', \n", |
| 189 | + " collections=collection_name,\n", |
| 190 | + " where=f\"exposure.day_obs = {day_obs} and exposure.seq_num >= {min_seq_num_30_sec} and exposure.seq_num <= {max_seq_num_30_sec} and instrument = 'LSSTComCam'\")" |
| 191 | + ] |
| 192 | + }, |
| 193 | + { |
| 194 | + "cell_type": "code", |
| 195 | + "execution_count": null, |
| 196 | + "id": "d32a7b53-e841-49d7-a6d0-5aa016e9f87e", |
| 197 | + "metadata": {}, |
| 198 | + "outputs": [], |
| 199 | + "source": [ |
| 200 | + "visit_dict_15_sec = dict()\n", |
| 201 | + "for visit_ref in visit_tables_15_sec:\n", |
| 202 | + " visit_table = butler.get(visit_ref)\n", |
| 203 | + " visit_dict_15_sec[visit_table.meta['visit']] = visit_table\n", |
| 204 | + "visit_dict_15_sec = dict(sorted(visit_dict_15_sec.items()))" |
| 205 | + ] |
| 206 | + }, |
| 207 | + { |
| 208 | + "cell_type": "code", |
| 209 | + "execution_count": null, |
| 210 | + "id": "ff407e43-7b1c-43b3-8ed5-f152380209e7", |
| 211 | + "metadata": {}, |
| 212 | + "outputs": [], |
| 213 | + "source": [ |
| 214 | + "visit_dict_30_sec = dict()\n", |
| 215 | + "for visit_ref in visit_tables_30_sec:\n", |
| 216 | + " visit_table = butler.get(visit_ref)\n", |
| 217 | + " visit_dict_30_sec[visit_table.meta['visit']] = visit_table\n", |
| 218 | + "visit_dict_30_sec = dict(sorted(visit_dict_30_sec.items()))" |
| 219 | + ] |
| 220 | + }, |
| 221 | + { |
| 222 | + "cell_type": "code", |
| 223 | + "execution_count": null, |
| 224 | + "id": "500b20a7-4371-4136-9063-66655a6636a9", |
| 225 | + "metadata": {}, |
| 226 | + "outputs": [], |
| 227 | + "source": [ |
| 228 | + "noll_idx = visit_table.meta['nollIndices']\n", |
| 229 | + "noll_min = np.min(noll_idx)\n", |
| 230 | + "noll_max = np.max(noll_idx)" |
| 231 | + ] |
| 232 | + }, |
| 233 | + { |
| 234 | + "cell_type": "code", |
| 235 | + "execution_count": null, |
| 236 | + "id": "848ade3f-9c12-47b1-96be-a218c4b107d7", |
| 237 | + "metadata": {}, |
| 238 | + "outputs": [], |
| 239 | + "source": [ |
| 240 | + "conv_array = getPsfGradPerZernike(jmin=noll_min, jmax=noll_max)" |
| 241 | + ] |
| 242 | + }, |
| 243 | + { |
| 244 | + "cell_type": "code", |
| 245 | + "execution_count": null, |
| 246 | + "id": "3cda025b-94fb-473d-bc73-d53a4cada57f", |
| 247 | + "metadata": {}, |
| 248 | + "outputs": [], |
| 249 | + "source": [ |
| 250 | + "zk_table_15_sec = Table()\n", |
| 251 | + "zk_table_15_sec = Table(names=['visit', 'full_array', 'fwhm_combined'], dtype=[int, (float, 25), float])\n", |
| 252 | + "for visit_id, visit_results in visit_dict_15_sec.items():\n", |
| 253 | + " zk_array = np.zeros((noll_max - noll_min + 1))\n", |
| 254 | + " for det_name in camera.getNameIter():\n", |
| 255 | + " zk_det_array = np.zeros((noll_max - noll_min + 1))\n", |
| 256 | + " zk_det_array[noll_idx - noll_min - 1] = visit_results[visit_results['detector'] == det_name]['zk_CCS']\n", |
| 257 | + " zk_array += zk_det_array * conv_array\n", |
| 258 | + " zk_array /= len(camera)\n", |
| 259 | + " zk_table_15_sec.add_row(vals=[visit_id, zk_array, np.sqrt(np.sum(zk_array**2))])" |
| 260 | + ] |
| 261 | + }, |
| 262 | + { |
| 263 | + "cell_type": "code", |
| 264 | + "execution_count": null, |
| 265 | + "id": "2d8f07c2-d467-4610-881b-bbb9c046b75b", |
| 266 | + "metadata": {}, |
| 267 | + "outputs": [], |
| 268 | + "source": [ |
| 269 | + "zk_table_30_sec = Table()\n", |
| 270 | + "zk_table_30_sec = Table(names=['visit', 'full_array', 'fwhm_combined'], dtype=[int, (float, 25), float])\n", |
| 271 | + "for visit_id, visit_results in visit_dict_30_sec.items():\n", |
| 272 | + " zk_array = np.zeros((noll_max - noll_min + 1))\n", |
| 273 | + " for det_name in camera.getNameIter():\n", |
| 274 | + " zk_det_array = np.zeros((noll_max - noll_min + 1))\n", |
| 275 | + " zk_det_array[noll_idx - noll_min - 1] = visit_results[visit_results['detector'] == det_name]['zk_CCS']\n", |
| 276 | + " zk_array += zk_det_array * conv_array\n", |
| 277 | + " zk_array /= len(camera)\n", |
| 278 | + " zk_table_30_sec.add_row(vals=[visit_id, zk_array, np.sqrt(np.sum(zk_array**2))])" |
| 279 | + ] |
| 280 | + }, |
| 281 | + { |
| 282 | + "cell_type": "code", |
| 283 | + "execution_count": null, |
| 284 | + "id": "c71d3a42-98f8-46ee-b2b5-13be0c569bf0", |
| 285 | + "metadata": {}, |
| 286 | + "outputs": [], |
| 287 | + "source": [ |
| 288 | + "plt.plot((zk_table_15_sec['visit'] - zk_table_15_sec['visit'][0])/3, zk_table_15_sec['fwhm_combined'], label='15 Second Exposures')\n", |
| 289 | + "plt.plot((zk_table_30_sec['visit'] - zk_table_30_sec['visit'][0])/3, zk_table_30_sec['fwhm_combined'], label='30 Second Exposures')\n", |
| 290 | + "plt.xlabel('Closed Loop Iteration')\n", |
| 291 | + "plt.ylabel('FWHM AOS Residual (arcsec)')\n", |
| 292 | + "plt.legend()\n", |
| 293 | + "band = visit_table.meta['band']\n", |
| 294 | + "plt.title(\n", |
| 295 | + " 'Closed Loop Convergence: 15 sec vs 30 sec.\\n ' +\n", |
| 296 | + " f'day_obs: {day_obs} seq_num: {min_seq_num_15_sec} - {max_seq_num_15_sec}, {min_seq_num_30_sec} - {max_seq_num_30_sec}, band: {band}'\n", |
| 297 | + ")" |
| 298 | + ] |
| 299 | + }, |
| 300 | + { |
| 301 | + "cell_type": "code", |
| 302 | + "execution_count": null, |
| 303 | + "id": "120b73b3-5cdf-42ff-a745-fe3ff8fe5e21", |
| 304 | + "metadata": {}, |
| 305 | + "outputs": [], |
| 306 | + "source": [] |
| 307 | + } |
| 308 | + ], |
| 309 | + "metadata": { |
| 310 | + "kernelspec": { |
| 311 | + "display_name": "LSST", |
| 312 | + "language": "python", |
| 313 | + "name": "lsst" |
| 314 | + }, |
| 315 | + "language_info": { |
| 316 | + "codemirror_mode": { |
| 317 | + "name": "ipython", |
| 318 | + "version": 3 |
| 319 | + }, |
| 320 | + "file_extension": ".py", |
| 321 | + "mimetype": "text/x-python", |
| 322 | + "name": "python", |
| 323 | + "nbconvert_exporter": "python", |
| 324 | + "pygments_lexer": "ipython3", |
| 325 | + "version": "3.12.9" |
| 326 | + } |
| 327 | + }, |
| 328 | + "nbformat": 4, |
| 329 | + "nbformat_minor": 5 |
| 330 | +} |
0 commit comments