Skip to content

Commit c1558e2

Browse files
authored
Merge pull request #1180 from dirac-institute/DES-Branch
DES branch
2 parents 5e9775c + 6429a94 commit c1558e2

40 files changed

+2996
-44
lines changed

docs/notebooks/DES_footprint.db

24 KB
Binary file not shown.
Lines changed: 297 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,297 @@
1+
{
2+
"cells": [
3+
{
4+
"cell_type": "markdown",
5+
"id": "a986e1d5",
6+
"metadata": {},
7+
"source": [
8+
"## Creating the SQLite Databases Necessary for DES Sorcha\n",
9+
"\n",
10+
"This notebook creates the SQLite databases needed to run **DES** in Sorcha. This code obtains the nescessory data from https://github.com/bernardinelli/DESTNOSIM/tree/master/data and converts it to sqlite databases. It creates two files:\n",
11+
"\n",
12+
"1. **Visits database** from `y6a1c.ccdcorners.fits.gz`\n",
13+
"2. **Pointing database** from `y6a1c.exposures.positions.fits`\n"
14+
]
15+
},
16+
{
17+
"cell_type": "code",
18+
"execution_count": null,
19+
"id": "70dc0f7f",
20+
"metadata": {},
21+
"outputs": [],
22+
"source": [
23+
"import numpy as np\n",
24+
"from astropy.io import fits\n",
25+
"from sqlalchemy import create_engine\n",
26+
"import sqlite3\n",
27+
"from astropy.time import Time\n",
28+
"import pandas as pd\n",
29+
"import requests\n",
30+
"from io import BytesIO\n",
31+
"import gzip\n",
32+
"\n",
33+
"\n",
34+
"\n",
35+
"def process_pointings_to_sqlite(\n",
36+
" use_url,\n",
37+
" fits_path=\"y6a1c.exposures.positions.fits\",\n",
38+
" db_path=\"DES_TNO.db\",\n",
39+
" url=\"https://github.com/bernardinelli/DESTNOSIM/raw/refs/heads/master/data/y6a1c.exposures.positions.fits\"\n",
40+
"):\n",
41+
" \"\"\"\n",
42+
" Processes a FITS file containing DES pointing data and stores \n",
43+
" data into a SQLite database with indexing.\n",
44+
"\n",
45+
" Parameters\n",
46+
" ----------\n",
47+
" use_url : bool\n",
48+
" If True, the FITS file will be downloaded from the specified `url`. If False, the local file path\n",
49+
" specified by `fits_path` will be used.\n",
50+
"\n",
51+
" fits_path : str, optional\n",
52+
" Path to the local FITS file (default is \"y6a1c.exposures.positions.fits\"). Only used if `use_url` is False.\n",
53+
"\n",
54+
" db_path : str, optional\n",
55+
" Path where the SQLite database will be saved (default is \"DES_TNO.db\").\n",
56+
"\n",
57+
" url : str, optional\n",
58+
" Direct URL to the raw FITS file hosted on GitHub.\n",
59+
"\n",
60+
" \n",
61+
" \"\"\"\n",
62+
" if use_url == True:\n",
63+
" response = requests.get(url)\n",
64+
" response.raise_for_status()\n",
65+
" with fits.open(BytesIO(response.content)) as hdul:\n",
66+
" df = hdul[1].data\n",
67+
" else:\n",
68+
" with fits.open(fits_path) as HDUl:\n",
69+
" df = HDUl[1].data\n",
70+
"\n",
71+
" # Define column names\n",
72+
" cov_1 = \"cov_xx\"\n",
73+
" cov_2 = \"cov_yy\"\n",
74+
" cov_3 = \"cov_xy\"\n",
75+
" observatory_1 = \"observatory_1\"\n",
76+
" observatory_2 = \"observatory_2\"\n",
77+
" observatory_3 = \"observatory_3\"\n",
78+
" velocity_1 = \"velocity_1\"\n",
79+
" velocity_2 = \"velocity_2\"\n",
80+
" velocity_3 = \"velocity_3\"\n",
81+
"\n",
82+
" # Create dictionary for DataFrame construction\n",
83+
" temp = {\n",
84+
" cov_1: np.array(df[\"cov\"][:, 0], dtype=\"<f8\"),\n",
85+
" cov_2: np.array(df[\"cov\"][:, 1], dtype=\"<f8\"),\n",
86+
" cov_3: np.array(df[\"cov\"][:, 2], dtype=\"<f8\"),\n",
87+
" \"covwarn\": df[\"covwarn\"],\n",
88+
" \"fieldDec\": np.array(df[\"dec\"], dtype=\"<f8\"),\n",
89+
" \"ecl_lat\": np.array(df[\"ecl_lat\"], dtype=\"<f8\"),\n",
90+
" \"ecl_lon\": np.array(df[\"ecl_lon\"], dtype=\"<f8\"),\n",
91+
" \"observationId\": np.array(df[\"expnum\"], \"<i4\"),\n",
92+
" \"filter\": df[\"filter\"],\n",
93+
" \"observationMidpointMJD\": np.array(df[\"mjd_mid\"], dtype=\"<f8\"),\n",
94+
" observatory_1: np.array(df[\"observatory\"][:, 0], dtype=\"<f8\"),\n",
95+
" observatory_2: np.array(df[\"observatory\"][:, 1], dtype=\"<f8\"),\n",
96+
" observatory_3: np.array(df[\"observatory\"][:, 2], dtype=\"<f8\"),\n",
97+
" \"obs_ecl_lon\": np.array(df[\"obs_ecl_lon\"], \"<f8\"),\n",
98+
" \"fieldRA\": np.array(df[\"ra\"], dtype=\"<f8\"),\n",
99+
" \"fiveSigmaDepth\": np.array(df[\"m50\"], dtype=\"<f8\"),\n",
100+
" \"k\": np.array(df[\"k\"], dtype=\"<f8\"),\n",
101+
" \"c\": np.array(df[\"c\"], dtype=\"<f8\"),\n",
102+
" velocity_1: np.array(df[\"velocity\"][:, 0], dtype=\"<f8\"),\n",
103+
" velocity_2: np.array(df[\"velocity\"][:, 1], dtype=\"<f8\"),\n",
104+
" velocity_3: np.array(df[\"velocity\"][:, 2], dtype=\"<f8\"),\n",
105+
" }\n",
106+
"\n",
107+
" # Set exposure times\n",
108+
" exo_time_s = np.full(len(df), 90)\n",
109+
" mask = (temp[\"observationMidpointMJD\"] < 57447) & (temp[\"filter\"] == \"Y\")\n",
110+
" exo_time_s[mask] = 45\n",
111+
"\n",
112+
" # Convert TDB to TAI\n",
113+
" time = Time(temp[\"observationMidpointMJD\"], format=\"mjd\", scale=\"utc\")\n",
114+
" time_TAI = time.tai\n",
115+
" temp[\"observationMidpointMJD\"] = time_TAI.value\n",
116+
"\n",
117+
" # Check for NaNs\n",
118+
" if np.any(pd.isnull(temp[\"observationMidpointMJD\"])):\n",
119+
" print(\"Warning: NaN values found in observationMidpointMJD\")\n",
120+
"\n",
121+
" # Create DataFrame and insert exposure times\n",
122+
" df_hdl1 = pd.DataFrame(temp)\n",
123+
" df_hdl1.insert(8, \"visitExposureTime\", exo_time_s)\n",
124+
"\n",
125+
" # Save to SQLite database\n",
126+
" engine = create_engine(f\"sqlite:///{db_path}\")\n",
127+
" df_hdl1.to_sql(\"observations\", engine, if_exists=\"replace\", index=False)\n",
128+
"\n",
129+
" # Create indexes efficiently\n",
130+
" index_queries = [\n",
131+
" \"CREATE INDEX idx_lat_long ON observations(ecl_lat,ecl_lon)\",\n",
132+
" \"CREATE INDEX idx_filter ON observations(filter)\",\n",
133+
" \"CREATE INDEX idx_dec_ra ON observations(fieldDec,fieldRA)\",\n",
134+
" \"CREATE INDEX idx_dec_ra_mjd ON observations(fieldDec,fieldRA,observationMidpointMJD)\",\n",
135+
" \"CREATE INDEX idx_mjd ON observations(observationMidpointMJD)\",\n",
136+
" \"CREATE INDEX idx_m50_c_k ON observations(fiveSigmaDepth,c,k)\"\n",
137+
" ]\n",
138+
"\n",
139+
" with sqlite3.connect(db_path) as db:\n",
140+
" cursor = db.cursor()\n",
141+
" for query in index_queries:\n",
142+
" cursor.execute(query)\n",
143+
"\n",
144+
" print(f\"DES pointing data processed and saved to: {db_path}\")\n",
145+
"\n",
146+
"\n",
147+
"\n",
148+
"def process_ccd_visits_to_sqlite(\n",
149+
" use_url,\n",
150+
" pointings_db_path,\n",
151+
" fits_path=\"y6a1c.ccdcorners.fits.gz\",\n",
152+
" db_path=\"DES_visits.db\",\n",
153+
" url=\"https://github.com/bernardinelli/DESTNOSIM/raw/refs/heads/master/data/y6a1c.ccdcorners.fits.gz\"\n",
154+
" ):\n",
155+
" \"\"\"\n",
156+
" Processes a FITS file containing DES ccd visits and stores \n",
157+
" data into a SQLite database with indexing.\n",
158+
"\n",
159+
" Parameters\n",
160+
" ----------\n",
161+
" use_url : bool\n",
162+
" If True, the FITS file will be downloaded from the specified `url`. If False, the local file path\n",
163+
" specified by `fits_path` will be used.\n",
164+
"\n",
165+
" fits_path : str, optional\n",
166+
" Path to the local FITS file (default is \"y6a1c.ccdcorners.fits.gz\"). \n",
167+
"\n",
168+
" db_path : str, optional\n",
169+
" Path where the SQLite database will be saved (default is \"DES_TNO.db\").\n",
170+
"\n",
171+
" url : str, optional\n",
172+
" Direct URL to the raw FITS file hosted on GitHub.\n",
173+
"\n",
174+
" \n",
175+
" \"\"\"\n",
176+
" if use_url == True:\n",
177+
" response = requests.get(url)\n",
178+
" response.raise_for_status()\n",
179+
" with gzip.open(BytesIO(response.content), 'rb') as gz:\n",
180+
" with fits.open(gz) as hdul:\n",
181+
" df = hdul[1].data\n",
182+
" else:\n",
183+
" with fits.open(fits_path) as HDUl:\n",
184+
" df = HDUl[1].data\n",
185+
"\n",
186+
"\n",
187+
" ra = df[\"ra\"]\n",
188+
" dec = df[\"dec\"]\n",
189+
"\n",
190+
" temp = {\n",
191+
" \"visitId\": df[\"expnum\"],\n",
192+
" \"detectorID\": df[\"ccdnum\"],\n",
193+
" \"llcra\": ra[:, 0],\n",
194+
" \"llcdec\": dec[:, 0],\n",
195+
" \"lrcra\": ra[:, 1],\n",
196+
" \"lrcdec\": dec[:, 1],\n",
197+
" \"urcra\": ra[:, 2],\n",
198+
" \"urcdec\": dec[:, 2],\n",
199+
" \"ulcra\": ra[:, 3],\n",
200+
" \"ulcdec\": dec[:, 3],\n",
201+
" \"ra\": ra[:, 4], # center RA\n",
202+
" \"dec\": dec[:, 4], # center Dec\n",
203+
" }\n",
204+
"\n",
205+
" # Create DataFrame\n",
206+
" df_hdl1 = pd.DataFrame(temp)\n",
207+
" \n",
208+
" engine_pointings = create_engine(f\"sqlite:///{pointings_db_path}\")\n",
209+
" df_depth = pd.read_sql(\"SELECT observationId, fiveSigmaDepth FROM observations\", engine_pointings)\n",
210+
" df_depth = df_depth.set_index(\"observationId\")\n",
211+
" df_depth = pd.read_sql(\"SELECT observationId, fiveSigmaDepth FROM observations\", engine_pointings)\n",
212+
" df_depth = df_depth.astype({\"fiveSigmaDepth\": \"<f8\"}) # just to be safe on dtype\n",
213+
" df_hdl1[\"magLim\"] = df_hdl1[\"visitId\"].map(df_depth[\"fiveSigmaDepth\"])\n",
214+
" \n",
215+
"\n",
216+
" # Save to SQLite\n",
217+
" engine = create_engine(f\"sqlite:///{db_path}\")\n",
218+
" df_hdl1.to_sql(\"observations\", engine, if_exists=\"replace\", index=False, chunksize=10000)\n",
219+
"\n",
220+
" # Create spatial index\n",
221+
" index_sql = \"\"\"\n",
222+
" CREATE INDEX IF NOT EXISTS idx_obs \n",
223+
" ON observations(visitId)\n",
224+
" \"\"\"\n",
225+
"\n",
226+
" with sqlite3.connect(db_path) as db:\n",
227+
" db.execute(index_sql)\n",
228+
"\n",
229+
" print(f\"DES CCD visits data processed and saved to {db_path}\")\n",
230+
"\n",
231+
"\n"
232+
]
233+
},
234+
{
235+
"cell_type": "code",
236+
"execution_count": 2,
237+
"id": "f09d3b0b",
238+
"metadata": {},
239+
"outputs": [
240+
{
241+
"name": "stdout",
242+
"output_type": "stream",
243+
"text": [
244+
"DES pointing data processed and saved to: DES_TNO.db\n",
245+
"DES CCD visits data processed and saved to DES_visits.db\n"
246+
]
247+
}
248+
],
249+
"source": [
250+
"pointing_db = \"DES_TNO.db\"\n",
251+
"process_pointings_to_sqlite(True,db_path=pointing_db)\n",
252+
"process_ccd_visits_to_sqlite(True,pointings_db_path=pointing_db)"
253+
]
254+
},
255+
{
256+
"cell_type": "markdown",
257+
"id": "ad8e5071",
258+
"metadata": {},
259+
"source": [
260+
"For running DES, obtain the databases using this notebook, run the command line arg `sorcha init` and select the **DES config file** option. Then run the code with:\n",
261+
"\n",
262+
"sorcha run -c DES_config_file.ini --pd DES_TNO.db --ob orbits_filename.csv -p colours_filename.csv -o ./ -s des --vd DES_visits.db\n",
263+
"\n",
264+
"`orbits_filename.csv` and `colours_filename.csv` should be your input files for object's orbits and colours, respectively."
265+
]
266+
},
267+
{
268+
"cell_type": "code",
269+
"execution_count": null,
270+
"id": "115468ba",
271+
"metadata": {},
272+
"outputs": [],
273+
"source": []
274+
}
275+
],
276+
"metadata": {
277+
"kernelspec": {
278+
"display_name": "sorcha",
279+
"language": "python",
280+
"name": "python3"
281+
},
282+
"language_info": {
283+
"codemirror_mode": {
284+
"name": "ipython",
285+
"version": 3
286+
},
287+
"file_extension": ".py",
288+
"mimetype": "text/x-python",
289+
"name": "python",
290+
"nbconvert_exporter": "python",
291+
"pygments_lexer": "ipython3",
292+
"version": "3.11.11"
293+
}
294+
},
295+
"nbformat": 4,
296+
"nbformat_minor": 5
297+
}

docs/notebooks/demo_VisitsFootprintFilter.ipynb

Lines changed: 320 additions & 0 deletions
Large diffs are not rendered by default.

docs/requirements.txt

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -21,3 +21,4 @@ sorcha-addons
2121
importlib_resources
2222
rubin-sim
2323
seaborn
24+
shapely

pyproject.toml

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -34,6 +34,7 @@ dependencies = [
3434
"pooch",
3535
"tqdm",
3636
"numba",
37+
"shapely",
3738
"importlib_resources",
3839
"scipy" # Needed for linear algebra in numba, do not remove!
3940
]

0 commit comments

Comments
 (0)