diff --git a/.gitignore b/.gitignore index 1f67b453..abb3a1ea 100644 --- a/.gitignore +++ b/.gitignore @@ -148,6 +148,8 @@ notebooks/de421.bsp notebooks/MPCORB.DAT.gz # Other temporary files +local/ +.clinerules/ tmp/* util/tmp/* .vscode/settings.json diff --git a/docs/conf.py b/docs/conf.py index 4ca84a36..71efa8e0 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -6,3 +6,4 @@ from documenteer.conf.guide import * # noqa: F403, import * linkcheck_retries = 2 +html_extra_path = ["extra"] diff --git a/docs/extra/visit_skymaps.html b/docs/extra/visit_skymaps.html new file mode 100644 index 00000000..ecab8d07 --- /dev/null +++ b/docs/extra/visit_skymaps.html @@ -0,0 +1,8922 @@ + + + + + +visit_skymaps + + + + + + + + + + + + +
+
+ +
+
+ +
+
+ +
+ +
+
+ +
+ +
+
+ +
+ +
+ + +
+
+ +
+ +
+
+ +
+ + +
+
+ +
+ +
+ + +
+
+ +
+ + +
+
+ +
+ + +
+
+ +
+ + +
+
+ +
+ + +
+
+ +
+ + +
+ +
+
+ + diff --git a/docs/tutorials.rst b/docs/tutorials.rst index 3a465c92..945da6d8 100644 --- a/docs/tutorials.rst +++ b/docs/tutorials.rst @@ -1,4 +1,4 @@ Tutorials ========= -For tutorial documentation, follow the tutorials in `schedview/notebooks/scheduler.ipynb` and `schedview/notebooks/prenight.ipynb`` in jupyter. +`Overview of the visit sky map builder `__. diff --git a/notebooks/architecture.ipynb b/notebooks/architecture.ipynb index c1925115..adca9d73 100644 --- a/notebooks/architecture.ipynb +++ b/notebooks/architecture.ipynb @@ -170,9 +170,7 @@ "import bokeh.models\n", "\n", "\n", - "def compute_minor_planet_positions(\n", - " minor_planet_orbits, start_mjd, end_mjd, time_step=7\n", - "):\n", + "def compute_minor_planet_positions(minor_planet_orbits, start_mjd, end_mjd, time_step=7):\n", " # Convert input fields into object appropriate for skyfield\n", " timescale = skyfield.api.load.timescale()\n", " start_ts = timescale.from_astropy(Time(start_mjd, format=\"mjd\"))\n", @@ -188,9 +186,7 @@ " orbit_rel_sun = skyfield.data.mpc.mpcorb_orbit(orbit, timescale, GM_SUN)\n", " minor_planet = sun + orbit_rel_sun\n", " for sample_time in sample_times:\n", - " ra, decl, distance = (\n", - " ephemeris[\"earth\"].at(sample_time).observe(minor_planet).radec()\n", - " )\n", + " ra, decl, distance = ephemeris[\"earth\"].at(sample_time).observe(minor_planet).radec()\n", " position_data[\"designation\"].append(orbit[\"designation\"])\n", " position_data[\"mjd\"].append(sample_time.to_astropy().mjd)\n", " position_data[\"ra\"].append(ra._degrees)\n", @@ -283,9 +279,7 @@ " factors=minor_planet_designations,\n", " )\n", "\n", - " figure.scatter(\n", - " \"ra\", \"decl\", color=cmap, legend_field=\"designation\", source=position_ds\n", - " )\n", + " figure.scatter(\"ra\", \"decl\", color=cmap, legend_field=\"designation\", source=position_ds)\n", " figure.title = \"Select minor planet positions\"\n", " figure.yaxis.axis_label = \"Declination (degrees)\"\n", " figure.xaxis.axis_label = \"R.A. (degrees)\"\n", @@ -394,9 +388,7 @@ " label=\"Start MJD\",\n", " )\n", "\n", - " end_mjd = param.Number(\n", - " default=60565, doc=\"Modified Julian Date of end of date window\", label=\"End MJD\"\n", - " )\n", + " end_mjd = param.Number(default=60565, doc=\"Modified Julian Date of end of date window\", label=\"End MJD\")\n", "\n", " orbits = param.Parameter()\n", "\n", diff --git a/notebooks/footprint_polygons.ipynb b/notebooks/footprint_polygons.ipynb index d9632418..e54bdc73 100644 --- a/notebooks/footprint_polygons.ipynb +++ b/notebooks/footprint_polygons.ipynb @@ -62,7 +62,7 @@ "metadata": {}, "outputs": [], "source": [ - "nside=64\n", + "nside = 64\n", "footprint_regions = get_current_footprint(nside)[1]\n", "footprint_polygons = schedview.compute.footprint.find_healpix_area_polygons(footprint_regions)\n", "footprint_polygons" @@ -101,7 +101,9 @@ "source": [ "psphere = Planisphere()\n", "asphere = ArmillarySphere()\n", - "schedview.plot.footprint.add_footprint_outlines_to_skymaps(footprint_polygons, [psphere, asphere], line_width=5)\n", + "schedview.plot.footprint.add_footprint_outlines_to_skymaps(\n", + " footprint_polygons, [psphere, asphere], line_width=5\n", + ")\n", "psphere.add_graticules(label_ra=True, label_decl=False)\n", "asphere.add_graticules()\n", "bokeh.io.show(bokeh.layouts.row([psphere.figure, asphere.figure]))" @@ -122,13 +124,13 @@ "metadata": {}, "outputs": [], "source": [ - "lowdust_polygons = footprint_polygons.loc['lowdust', :]\n", + "lowdust_polygons = footprint_polygons.loc[\"lowdust\", :]\n", "psphere = Planisphere()\n", "asphere = ArmillarySphere()\n", "schedview.plot.footprint.add_footprint_outlines_to_skymaps(lowdust_polygons, [psphere, asphere], line_width=5)\n", "psphere.add_graticules(label_ra=True, label_decl=False)\n", "asphere.add_graticules()\n", - "bokeh.io.show(bokeh.layouts.row([psphere.figure, asphere.figure]))\n" + "bokeh.io.show(bokeh.layouts.row([psphere.figure, asphere.figure]))" ] }, { @@ -146,11 +148,15 @@ "metadata": {}, "outputs": [], "source": [ - "simplified_footprint_polygons = schedview.compute.footprint.find_healpix_area_polygons(footprint_regions, simplify_tolerance=0.01)\n", + "simplified_footprint_polygons = schedview.compute.footprint.find_healpix_area_polygons(\n", + " footprint_regions, simplify_tolerance=0.01\n", + ")\n", "print(f\"Vertexes used dropped from {len(footprint_polygons)} to {len(simplified_footprint_polygons)}\")\n", "psphere = Planisphere()\n", "asphere = ArmillarySphere()\n", - "schedview.plot.footprint.add_footprint_outlines_to_skymaps(simplified_footprint_polygons, [psphere, asphere], line_width=5)\n", + "schedview.plot.footprint.add_footprint_outlines_to_skymaps(\n", + " simplified_footprint_polygons, [psphere, asphere], line_width=5\n", + ")\n", "psphere.add_graticules(label_ra=True, label_decl=False)\n", "asphere.add_graticules()\n", "bokeh.io.show(bokeh.layouts.row([psphere.figure, asphere.figure]))" @@ -174,12 +180,14 @@ "outputs": [], "source": [ "psphere = Planisphere()\n", - "psphere.plot.title=\"Lambert Azimuthal Equal Area projection\"\n", + "psphere.plot.title = \"Lambert Azimuthal Equal Area projection\"\n", "asphere = ArmillarySphere()\n", - "asphere.plot.title=\"Dynamic Orthographic projection\"\n", + "asphere.plot.title = \"Dynamic Orthographic projection\"\n", "msphere = MollweideMap()\n", - "msphere.plot.title=\"Mollweide projection\"\n", - "schedview.plot.footprint.add_footprint_outlines_to_skymaps(footprint_polygons, [psphere, asphere, msphere], filled=True)\n", + "msphere.plot.title = \"Mollweide projection\"\n", + "schedview.plot.footprint.add_footprint_outlines_to_skymaps(\n", + " footprint_polygons, [psphere, asphere, msphere], filled=True\n", + ")\n", "psphere.add_graticules(label_ra=True, label_decl=False)\n", "asphere.add_graticules()\n", "msphere.add_graticules()\n", diff --git a/notebooks/introduction.ipynb b/notebooks/introduction.ipynb index acc8beda..e000ce52 100644 --- a/notebooks/introduction.ipynb +++ b/notebooks/introduction.ipynb @@ -512,7 +512,9 @@ "\n", "from rubin_sim import maf\n", "\n", - "visits_fname = \"/sdf/group/rubin/web_data/sim-data/sims_featureScheduler_runs3.5/baseline/baseline_v3.5_10yrs.db\"\n", + "visits_fname = (\n", + " \"/sdf/group/rubin/web_data/sim-data/sims_featureScheduler_runs3.5/baseline/baseline_v3.5_10yrs.db\"\n", + ")\n", "\n", "# In practice you'll need to convert day_obs to opsim_night, not shown here.\n", "opsim_night = 10\n", @@ -655,13 +657,9 @@ "os.environ[\"LSST_DISABLE_BUCKET_VALIDATION\"] = \"1\"\n", "os.environ[\"S3_ENDPOINT_URL\"] = \"https://s3dfrgw.slac.stanford.edu/\"\n", "sim_archive_rp = (\n", - " ResourcePath(archive_uri)\n", - " .join(sim_date, forceDirectory=True)\n", - " .join(f\"{sim_index}\", forceDirectory=True)\n", - ")\n", - "sim_archive_metadata = yaml.safe_load(\n", - " sim_archive_rp.join(\"sim_metadata.yaml\").read().decode()\n", + " ResourcePath(archive_uri).join(sim_date, forceDirectory=True).join(f\"{sim_index}\", forceDirectory=True)\n", ")\n", + "sim_archive_metadata = yaml.safe_load(sim_archive_rp.join(\"sim_metadata.yaml\").read().decode())\n", "sim_rp = sim_archive_rp.join(sim_archive_metadata[\"files\"][\"observations\"][\"name\"])\n", "sim_rp" ] @@ -861,7 +859,9 @@ "metadata": {}, "outputs": [], "source": [ - "visit_display = f'
{visits.to_html()}
'\n", + "visit_display = (\n", + " f'
{visits.to_html()}
'\n", + ")\n", "display(HTML(visit_display))" ] }, diff --git a/notebooks/make_test_schedulers.ipynb b/notebooks/make_test_schedulers.ipynb index c946603f..c4002932 100644 --- a/notebooks/make_test_schedulers.ipynb +++ b/notebooks/make_test_schedulers.ipynb @@ -69,7 +69,7 @@ }, "outputs": [], "source": [ - "SchedulerPickleContent = namedtuple('SchedulerPickleContent', ['scheduler', 'conditions'])" + "SchedulerPickleContent = namedtuple(\"SchedulerPickleContent\", [\"scheduler\", \"conditions\"])" ] }, { @@ -84,7 +84,7 @@ "mjd_start = survey_start_mjd()\n", "model_observatory = ModelObservatory(mjd_start=mjd_start)\n", "nside = model_observatory.nside\n", - "sim_duration = (4*u.hour).to(u.day).value" + "sim_duration = (4 * u.hour).to(u.day).value" ] }, { @@ -186,21 +186,16 @@ }, "outputs": [], "source": [ - "def make_sky_bf_list(band='g', nside=32):\n", + "def make_sky_bf_list(band=\"g\", nside=32):\n", " not_twilight = bf.feasibility_funcs.NotTwilightBasisFunction()\n", " moon_limit = bf.mask_basis_funcs.MoonAvoidanceBasisFunction(nside=nside, moon_distance=30.0)\n", " zenith_limit = bf.mask_basis_funcs.ZenithShadowMaskBasisFunction(nside=nside, min_alt=20.0, max_alt=82.0)\n", - " sky_brightness_limit = bf.basis_functions.SkybrightnessLimitBasisFunction(nside=nside, filtername=band, sbmin=18.5, sbmax=30)\n", + " sky_brightness_limit = bf.basis_functions.SkybrightnessLimitBasisFunction(\n", + " nside=nside, filtername=band, sbmin=18.5, sbmax=30\n", + " )\n", " wind_limit = bf.basis_functions.AvoidDirectWind(5)\n", " m5diff = bf.basis_functions.M5DiffBasisFunction(filtername=band, nside=nside)\n", - " basis_functions = [\n", - " not_twilight,\n", - " moon_limit,\n", - " zenith_limit,\n", - " sky_brightness_limit,\n", - " wind_limit,\n", - " m5diff\n", - " ]\n", + " basis_functions = [not_twilight, moon_limit, zenith_limit, sky_brightness_limit, wind_limit, m5diff]\n", " return basis_functions" ] }, @@ -213,9 +208,11 @@ }, "outputs": [], "source": [ - "def make_field_bf_list(ra, decl, band='g', nside=32):\n", + "def make_field_bf_list(ra, decl, band=\"g\", nside=32):\n", " basis_functions = make_sky_bf_list(band=band, nside=nside)\n", - " basis_functions.append(bf.feasibility_funcs.HourAngleLimitBasisFunction(RA=ra, ha_limits=[[22,24], [0,2]]))\n", + " basis_functions.append(\n", + " bf.feasibility_funcs.HourAngleLimitBasisFunction(RA=ra, ha_limits=[[22, 24], [0, 2]])\n", + " )\n", " return basis_functions" ] }, @@ -228,12 +225,21 @@ }, "outputs": [], "source": [ - "def make_field_survey(ra, decl, band='g', nside=32):\n", + "def make_field_survey(ra, decl, band=\"g\", nside=32):\n", " basis_functions = make_field_bf_list(ra, decl, band=band, nside=nside)\n", " sequence = band\n", - " nvis = [1]*len(band)\n", + " nvis = [1] * len(band)\n", " survey_name = f\"field_{ra}_{'n' if decl<0 else 'p'}{np.abs(decl)}_{band}\"\n", - " survey = FieldSurvey(basis_functions, np.array([ra]), np.array([decl]), sequence=sequence, nvis=nvis, nside=nside, survey_name=survey_name, reward_value=1.0)\n", + " survey = FieldSurvey(\n", + " basis_functions,\n", + " np.array([ra]),\n", + " np.array([decl]),\n", + " sequence=sequence,\n", + " nvis=nvis,\n", + " nside=nside,\n", + " survey_name=survey_name,\n", + " reward_value=1.0,\n", + " )\n", " return survey" ] }, @@ -263,7 +269,7 @@ "outputs": [], "source": [ "decl = 0\n", - "band = 'g'\n", + "band = \"g\"\n", "field_surveys = [make_field_survey(ra, decl, band, nside) for ra in range(0, 360, 15)]" ] }, @@ -286,7 +292,11 @@ "source": [ "sky_basis_functions = make_sky_bf_list(band=band, nside=nside)\n", "weights = [1] * len(sky_basis_functions)\n", - "greedy_surveys = [rubin_scheduler.scheduler.surveys.surveys.GreedySurvey(sky_basis_functions, weights, filtername=band, survey_name=f\"greedy_{band}\")]" + "greedy_surveys = [\n", + " rubin_scheduler.scheduler.surveys.surveys.GreedySurvey(\n", + " sky_basis_functions, weights, filtername=band, survey_name=f\"greedy_{band}\"\n", + " )\n", + "]" ] }, { @@ -374,10 +384,10 @@ "outputs": [], "source": [ "output_tuple = (scheduler, conditions)\n", - "if not os.path.isdir('../tmp'):\n", - " os.mkdir('../tmp')\n", - "fname = Path('../tmp/eq_field_survey_v0.p.xz').resolve()\n", - "with lzma.open(fname, 'wb') as sched_out:\n", + "if not os.path.isdir(\"../tmp\"):\n", + " os.mkdir(\"../tmp\")\n", + "fname = Path(\"../tmp/eq_field_survey_v0.p.xz\").resolve()\n", + "with lzma.open(fname, \"wb\") as sched_out:\n", " pickle.dump(output_tuple, sched_out)\n", "\n", "fname" @@ -413,7 +423,7 @@ "outputs": [], "source": [ "scheduler_reward_df = scheduler.make_reward_df(conditions, accum=True)\n", - "scheduler_reward_df.loc[(0,1),:]" + "scheduler_reward_df.loc[(0, 1), :]" ] }, { @@ -472,10 +482,12 @@ }, "outputs": [], "source": [ - "month_mjds = np.arange(int(mjd_start), int(mjd_start)+30)\n", - "month_phases = np.array([model_observatory.almanac.get_sun_moon_positions(mjd)['moon_phase'] for mjd in month_mjds])\n", + "month_mjds = np.arange(int(mjd_start), int(mjd_start) + 30)\n", + "month_phases = np.array(\n", + " [model_observatory.almanac.get_sun_moon_positions(mjd)[\"moon_phase\"] for mjd in month_mjds]\n", + ")\n", "mjd = month_mjds[np.argmax(month_phases)]\n", - "Time(mjd, format='mjd').iso" + "Time(mjd, format=\"mjd\").iso" ] }, { @@ -497,12 +509,12 @@ "source": [ "sunset_info = model_observatory.almanac.get_sunset_info(mjd)\n", "sunset, sunrise = sunset_info[3], sunset_info[4]\n", - "for mjd in np.arange(sunset, sunrise, 1.0/(24*4)):\n", - " moon_alt = model_observatory.almanac.get_sun_moon_positions(mjd)['moon_alt']\n", - " if np.abs(np.degrees(moon_alt)-55) < 5:\n", + "for mjd in np.arange(sunset, sunrise, 1.0 / (24 * 4)):\n", + " moon_alt = model_observatory.almanac.get_sun_moon_positions(mjd)[\"moon_alt\"]\n", + " if np.abs(np.degrees(moon_alt) - 55) < 5:\n", " break\n", "\n", - "time_to_sched = Time(mjd, format='mjd')\n", + "time_to_sched = Time(mjd, format=\"mjd\")\n", "time_to_sched.iso, mjd" ] }, @@ -535,8 +547,9 @@ }, "outputs": [], "source": [ - "moon_az = model_observatory.almanac.get_sun_moon_positions(time_to_sched.mjd)['moon_az']\n", - "wind_az = moon_az - np.pi ;# opposite moon\n", + "moon_az = model_observatory.almanac.get_sun_moon_positions(time_to_sched.mjd)[\"moon_az\"]\n", + "wind_az = moon_az - np.pi\n", + "# opposite moon\n", "wind_data = rubin_scheduler.site_models.ConstantWindData(wind_speed=18.0, wind_direction=wind_az)\n", "model_observatory.wind_data = wind_data" ] @@ -593,8 +606,8 @@ "outputs": [], "source": [ "output_tuple = (scheduler, conditions)\n", - "fname = Path('../tmp/eq_field_survey_v0_infeasible1.p.xz').resolve()\n", - "with lzma.open(fname, 'wb') as sched_out:\n", + "fname = Path(\"../tmp/eq_field_survey_v0_infeasible1.p.xz\").resolve()\n", + "with lzma.open(fname, \"wb\") as sched_out:\n", " pickle.dump(output_tuple, sched_out)\n", "\n", "fname" diff --git a/notebooks/prenight.ipynb b/notebooks/prenight.ipynb index 2b445dbe..296c55f2 100644 --- a/notebooks/prenight.ipynb +++ b/notebooks/prenight.ipynb @@ -297,10 +297,7 @@ "outputs": [], "source": [ "# If the date represents the local calendar date at sunset, we need to shift by the longitude in units of days\n", - "this_night = (\n", - " np.floor(observatory.almanac.sunsets[\"sunset\"] + observatory.site.longitude / 360)\n", - " == evening_mjd\n", - ")\n", + "this_night = np.floor(observatory.almanac.sunsets[\"sunset\"] + observatory.site.longitude / 360) == evening_mjd\n", "\n", "sim_start_mjd = observatory.almanac.sunsets[this_night][\"sun_n12_setting\"][0]\n", "sim_end_mjd = observatory.almanac.sunsets[this_night][\"sunrise\"][0]\n", @@ -423,9 +420,7 @@ "metadata": {}, "outputs": [], "source": [ - "with NamedTemporaryFile(\n", - " prefix=\"scheduler-\", suffix=\".pickle.xz\", dir=data_dir.name\n", - ") as temp_file:\n", + "with NamedTemporaryFile(prefix=\"scheduler-\", suffix=\".pickle.xz\", dir=data_dir.name) as temp_file:\n", " scheduler_fname = temp_file.name\n", "\n", "with lzma.open(scheduler_fname, \"wb\", format=lzma.FORMAT_XZ) as pio:\n", @@ -441,9 +436,7 @@ "metadata": {}, "outputs": [], "source": [ - "with NamedTemporaryFile(\n", - " prefix=\"rewards-\", suffix=\".h5\", dir=data_dir.name\n", - ") as temp_file:\n", + "with NamedTemporaryFile(prefix=\"rewards-\", suffix=\".h5\", dir=data_dir.name) as temp_file:\n", " rewards_fname = temp_file.name\n", "\n", "reward_df.to_hdf(rewards_fname, \"reward_df\")\n", @@ -567,9 +560,7 @@ }, "outputs": [], "source": [ - "with NamedTemporaryFile(\n", - " prefix=\"custom_prenight_tabs-\", suffix=\".json\", dir=data_dir.name\n", - ") as temp_file:\n", + "with NamedTemporaryFile(prefix=\"custom_prenight_tabs-\", suffix=\".json\", dir=data_dir.name) as temp_file:\n", " custom_tabs_fname = temp_file.name\n", "\n", "with open(custom_tabs_fname, \"w\") as custom_tabs_file:\n", diff --git a/notebooks/prenight_matplotlib_extension.ipynb b/notebooks/prenight_matplotlib_extension.ipynb index e5df3fa5..31c77943 100644 --- a/notebooks/prenight_matplotlib_extension.ipynb +++ b/notebooks/prenight_matplotlib_extension.ipynb @@ -235,9 +235,7 @@ "\n", " # Add your new plot to this dictionary\n", " new_tab_name = \"Azimuth wrap\"\n", - " tab_contents[new_tab_name] = pn.param.ParamMethod(\n", - " self.make_az_wrap_plot, loading_indicator=True\n", - " )\n", + " tab_contents[new_tab_name] = pn.param.ParamMethod(self.make_az_wrap_plot, loading_indicator=True)\n", "\n", " return tab_contents" ] @@ -294,11 +292,7 @@ " opsim = schedview.collect.opsim.read_opsim(opsim_fname)\n", " start_mjd = opsim.observationStartMJD.min()\n", " start_datetime_utc = Time(start_mjd, format=\"mjd\").datetime\n", - " night_date = (\n", - " pd.Timestamp(start_datetime_utc, tz=\"UTC\")\n", - " .tz_convert(\"Chile/Continental\")\n", - " .date()\n", - " )\n", + " night_date = pd.Timestamp(start_datetime_utc, tz=\"UTC\").tz_convert(\"Chile/Continental\").date()\n", " return night_date\n", "\n", "\n", @@ -614,9 +608,7 @@ " if self._reward_df is None or self._obs_rewards is None:\n", " return \"No rewards are loaded\"\n", "\n", - " return pn.pane.Matplotlib(\n", - " my_plot(self._visits, self._reward_df, self._obs_rewards)\n", - " )\n", + " return pn.pane.Matplotlib(my_plot(self._visits, self._reward_df, self._obs_rewards))\n", "\n", " def initialize_tab_contents(self):\n", " # Start with the dictionary with the tabs defined\n", @@ -625,9 +617,7 @@ "\n", " # Add your new plot to this dictionary\n", " new_tab_name = \"My plot\"\n", - " tab_contents[new_tab_name] = pn.param.ParamMethod(\n", - " self.make_my_plot, loading_indicator=True\n", - " )\n", + " tab_contents[new_tab_name] = pn.param.ParamMethod(self.make_my_plot, loading_indicator=True)\n", "\n", " return tab_contents" ] diff --git a/notebooks/prenight_multielement_extension.ipynb b/notebooks/prenight_multielement_extension.ipynb index f7e46e47..656cd548 100644 --- a/notebooks/prenight_multielement_extension.ipynb +++ b/notebooks/prenight_multielement_extension.ipynb @@ -218,9 +218,7 @@ " \"sunAlt\": \"max sun alt\",\n", " }\n", "\n", - " table_widget = pn.widgets.Tabulator(\n", - " activity_summary, titles=titles, show_index=False\n", - " )\n", + " table_widget = pn.widgets.Tabulator(activity_summary, titles=titles, show_index=False)\n", " return table_widget" ] }, @@ -305,9 +303,7 @@ " new_tab_name = \"Activities\"\n", " tab_contents[new_tab_name] = pn.Row(\n", " pn.param.ParamMethod(self.make_activity_bars, loading_indicator=True),\n", - " pn.param.ParamMethod(\n", - " self.make_activity_summary_table, loading_indicator=True\n", - " ),\n", + " pn.param.ParamMethod(self.make_activity_summary_table, loading_indicator=True),\n", " )\n", "\n", " return tab_contents" @@ -364,11 +360,7 @@ " opsim = schedview.collect.opsim.read_opsim(opsim_fname)\n", " start_mjd = opsim.observationStartMJD.min()\n", " start_datetime_utc = Time(start_mjd, format=\"mjd\").datetime\n", - " night_date = (\n", - " pd.Timestamp(start_datetime_utc, tz=\"UTC\")\n", - " .tz_convert(\"Chile/Continental\")\n", - " .date()\n", - " )\n", + " night_date = pd.Timestamp(start_datetime_utc, tz=\"UTC\").tz_convert(\"Chile/Continental\").date()\n", " return night_date\n", "\n", "\n", diff --git a/notebooks/prenight_rewardplot.ipynb b/notebooks/prenight_rewardplot.ipynb index ed4589e3..dfd67499 100644 --- a/notebooks/prenight_rewardplot.ipynb +++ b/notebooks/prenight_rewardplot.ipynb @@ -245,9 +245,7 @@ "outputs": [], "source": [ "def configure_survey_selector(survey_selector, reward_df, tier):\n", - " surveys = (\n", - " reward_df.set_index(\"tier_label\").loc[tier, \"survey_label\"].unique().tolist()\n", - " )\n", + " surveys = reward_df.set_index(\"tier_label\").loc[tier, \"survey_label\"].unique().tolist()\n", " survey_selector.options = surveys\n", " survey_selector.value = surveys[:10] if len(surveys) > 10 else surveys" ] @@ -262,9 +260,7 @@ "outputs": [], "source": [ "def configure_basis_function_selector(basis_function_selector, reward_df, tier):\n", - " basis_functions = (\n", - " reward_df.set_index(\"tier_label\").loc[tier, \"basis_function\"].unique().tolist()\n", - " )\n", + " basis_functions = reward_df.set_index(\"tier_label\").loc[tier, \"basis_function\"].unique().tolist()\n", " basis_function_selector.options = [\"Total\"] + basis_functions\n", " basis_function_selector.value = \"Total\"" ] @@ -286,9 +282,7 @@ " width_policy=\"fit\",\n", ")\n", "\n", - "surveys = pn.widgets.MultiSelect(\n", - " name=\"Displayed surveys\", options=[\"foo\"], value=[\"foo\"], width_policy=\"fit\"\n", - ")\n", + "surveys = pn.widgets.MultiSelect(name=\"Displayed surveys\", options=[\"foo\"], value=[\"foo\"], width_policy=\"fit\")\n", "configure_survey_selector(surveys, reward_df, tier.value)\n", "\n", "\n", diff --git a/notebooks/scheduler.ipynb b/notebooks/scheduler.ipynb index dae9d27a..74f733f4 100644 --- a/notebooks/scheduler.ipynb +++ b/notebooks/scheduler.ipynb @@ -280,7 +280,7 @@ "\n", "if remake_scheduler_pickle:\n", " # Use standard start time so standard skybrightness_pre small file works\n", - " start_date = Time(survey_start_mjd(), format='mjd', scale='utc') + TimeDelta(1, format='jd')\n", + " start_date = Time(survey_start_mjd(), format=\"mjd\", scale=\"utc\") + TimeDelta(1, format=\"jd\")\n", " start_date = start_date.isot.split(\"T\")[0]\n", " print(f\"Using start date of {start_date}\")\n", " local_timezone = ZoneInfo(\"Chile/Continental\")\n", @@ -312,9 +312,7 @@ }, "outputs": [], "source": [ - "scheduler, conditions = schedview.collect.scheduler_pickle.read_scheduler(\n", - " scheduler_pickle\n", - ")" + "scheduler, conditions = schedview.collect.scheduler_pickle.read_scheduler(scheduler_pickle)" ] }, { @@ -376,7 +374,7 @@ }, "outputs": [], "source": [ - "new_time = start_time + TimeDelta(0.25, format='jd')\n", + "new_time = start_time + TimeDelta(0.25, format=\"jd\")\n", "# Just check we didn't accidentally go past sunrise\n", "assert new_time.mjd < conditions.sun_n18_rising" ] @@ -676,9 +674,7 @@ }, "outputs": [], "source": [ - "schedview.compute.survey.make_survey_reward_df(\n", - " survey, conditions, reward_df.loc[[(tier_id, survey_id)], :]\n", - ")" + "schedview.compute.survey.make_survey_reward_df(survey, conditions, reward_df.loc[[(tier_id, survey_id)], :])" ] }, { @@ -753,9 +749,7 @@ "outputs": [], "source": [ "healpix_key = \"g_sky\"\n", - "sky_map = schedview.plot.survey.map_survey_healpix(\n", - " conditions.mjd, survey_maps, healpix_key, nside\n", - ")\n", + "sky_map = schedview.plot.survey.map_survey_healpix(conditions.mjd, survey_maps, healpix_key, nside)\n", "sky_map.notebook_display()" ] }, @@ -785,12 +779,8 @@ }, "outputs": [], "source": [ - "def show_scheduler(\n", - " pickle_fname, mjd, tier_id=2, survey_id=2, map_key=\"reward\", nside=None\n", - "):\n", - " scheduler, conditions = schedview.collect.scheduler_pickle.read_scheduler(\n", - " pickle_fname\n", - " )\n", + "def show_scheduler(pickle_fname, mjd, tier_id=2, survey_id=2, map_key=\"reward\", nside=None):\n", + " scheduler, conditions = schedview.collect.scheduler_pickle.read_scheduler(pickle_fname)\n", "\n", " if nside is None:\n", " nside = conditions.nside\n", @@ -806,9 +796,7 @@ " survey_maps = schedview.compute.survey.compute_maps(survey, conditions, nside=nside)\n", "\n", " # Display the selected map\n", - " sky_map = schedview.plot.survey.map_survey_healpix(\n", - " conditions.mjd, survey_maps, healpix_key, nside\n", - " )\n", + " sky_map = schedview.plot.survey.map_survey_healpix(conditions.mjd, survey_maps, healpix_key, nside)\n", " sky_map.notebook_display()\n", "\n", " # Display the selected survey basis function summary\n", @@ -824,9 +812,7 @@ " display(scheduler_summary_df)\n", "\n", "\n", - "mjd = Time(\n", - " pd.Timestamp(f\"{start_date} 23:00:00\", tzinfo=ZoneInfo(\"Chile/Continental\"))\n", - ").mjd\n", + "mjd = Time(pd.Timestamp(f\"{start_date} 23:00:00\", tzinfo=ZoneInfo(\"Chile/Continental\"))).mjd\n", "show_scheduler(scheduler_pickle, mjd, 2, 2, \"reward\")" ] }, diff --git a/notebooks/visit_skymaps.ipynb b/notebooks/visit_skymaps.ipynb new file mode 100644 index 00000000..ed76c134 --- /dev/null +++ b/notebooks/visit_skymaps.ipynb @@ -0,0 +1,577 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "0", + "metadata": {}, + "source": [ + "# About `schedview.plot.visit_skymaps`\n", + "\n", + "This notebook demonstrates the use of `schedview.plot.visit_skymaps.VisitMapBuilder`.\n", + "\n", + "To convert this page to an html file: `jupyter nbconvert --to html --execute visit_skymaps.ipynb`" + ] + }, + { + "cell_type": "markdown", + "id": "1", + "metadata": {}, + "source": [ + "## Imports and notebook setup" + ] + }, + { + "cell_type": "markdown", + "id": "2", + "metadata": {}, + "source": [ + "### Imports of general infrastructure" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "3", + "metadata": {}, + "outputs": [], + "source": [ + "import sys\n", + "import os\n", + "from pathlib import Path\n", + "\n", + "import bokeh.io\n", + "import numpy as np" + ] + }, + { + "cell_type": "markdown", + "id": "4", + "metadata": {}, + "source": [ + "### Adjust `PYTHONPATH` to point to development directories for `schedview` and `uranography`, if desired\n", + "\n", + "If you want to use the version provided by the current kernel, just leave `devel_dir` as `None`.\n", + "Otherwise, update `devel_dir` as appropriate." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "5", + "metadata": {}, + "outputs": [], + "source": [ + "devel_dir = None\n", + "# devel_dir = Path('/sdf/data/rubin/user/neilsen/devel')\n", + "\n", + "if devel_dir is not None:\n", + " for module_name in (\"schedview\", \"uranography\", \"rubin_sim\"):\n", + " module_path = devel_dir.joinpath(module_name)\n", + " if module_path.exists():\n", + " sys.path.insert(0, module_path.as_posix())" + ] + }, + { + "cell_type": "markdown", + "id": "6", + "metadata": {}, + "source": [ + "### Import project specific modules\n", + "\n", + "Use `autoreload` for `uranography` and `schedview` to support development." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "7", + "metadata": {}, + "outputs": [], + "source": [ + "%load_ext autoreload\n", + "%autoreload 1" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "8", + "metadata": {}, + "outputs": [], + "source": [ + "from rubin_sim.data import get_baseline\n", + "from rubin_scheduler.scheduler.utils import get_current_footprint\n", + " \n", + "%aimport uranography.spheremap\n", + "%aimport uranography.planisphere\n", + "%aimport uranography.armillary\n", + "%aimport uranography.camera\n", + "%aimport uranography.horizon\n", + "%aimport uranography.mollweide\n", + "%aimport schedview.plot.visit_skymaps\n", + "%aimport schedview.collect" + ] + }, + { + "cell_type": "markdown", + "id": "9", + "metadata": {}, + "source": [ + "### Set basic parameters" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "10", + "metadata": {}, + "outputs": [], + "source": [ + "# Night of the baseline to use for example data.\n", + "# For example, 365 is one year into the survey.\n", + "night = 365\n", + "\n", + "# healpy nside for maps to be shown.\n", + "nside = 64" + ] + }, + { + "cell_type": "markdown", + "id": "11", + "metadata": {}, + "source": [ + "### Prepare notebook for bokeh output" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "12", + "metadata": {}, + "outputs": [], + "source": [ + "bokeh.io.output_notebook()" + ] + }, + { + "cell_type": "markdown", + "id": "13", + "metadata": {}, + "source": [ + "### Read sample data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "14", + "metadata": {}, + "outputs": [], + "source": [ + "visits = (\n", + " schedview.collect.read_opsim(get_baseline())\n", + " .reset_index()\n", + " .set_index(\"night\")\n", + " .loc[[night, night + 1, night + 2], :]\n", + " .reset_index()\n", + " .set_index(\"observationId\")\n", + " .copy()\n", + ")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "15", + "metadata": {}, + "outputs": [], + "source": [ + "footprint_depth_by_band, footprint_regions = get_current_footprint(nside)" + ] + }, + { + "cell_type": "markdown", + "id": "16", + "metadata": {}, + "source": [ + "## Minimal example\n", + "\n", + "The `VisitMapBuilder` uses a \"[fluent method chaining](https://en.wikipedia.org/wiki/Fluent_interface)\" [builder pattern](https://en.wikipedia.org/wiki/Builder_pattern) interface: the calling function creates an instance of `VisitMapBuilder`, followed by a sequence of method calls that set the various desired features and parameters, and finally a call to `build` which returns an instance of the desired object. In this case, the object returned by `build` is an instance of `bokeh.models.UIElement`, a `bokeh` object that can be rendered with `bokeh.io` functions to a file, in an application interface or on a web page, or in a jupyter notebook." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "17", + "metadata": {}, + "outputs": [], + "source": [ + "builder = (\n", + " schedview.plot.visit_skymaps.VisitMapBuilder(\n", + " map_classes=[uranography.armillary.ArmillarySphere, uranography.planisphere.Planisphere],\n", + " )\n", + " .add_graticules()\n", + ")\n", + "viewable = builder.build()\n", + "bokeh.io.show(viewable)" + ] + }, + { + "cell_type": "markdown", + "id": "18", + "metadata": {}, + "source": [ + "## Example with many elements\n", + "\n", + "Methods can be called to provide additional data, features and controls:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "19", + "metadata": {}, + "outputs": [], + "source": [ + "builder = (\n", + " schedview.plot.visit_skymaps.VisitMapBuilder(\n", + " visits,\n", + " mjd=visits[\"observationStartMJD\"].max(),\n", + " map_classes=[uranography.armillary.ArmillarySphere, uranography.planisphere.Planisphere],\n", + " figure_kwargs={\"match_aspect\": True},\n", + " )\n", + " .add_visit_patches()\n", + " .add_footprint_outlines(footprint_regions)\n", + " .hide_horizon_sliders()\n", + " .make_up_north()\n", + " .add_eq_sliders()\n", + " .add_graticules()\n", + " .add_ecliptic()\n", + " .add_galactic_plane()\n", + " .add_datetime_slider()\n", + " .hide_mjd_slider()\n", + " .hide_future_and_other_night_visits()\n", + " .highlight_recent_visits()\n", + " .add_body(\"sun\", size=15, color=\"yellow\", alpha=1.0)\n", + " .add_body(\"moon\", size=15, color=\"orange\", alpha=0.8)\n", + " .add_horizon()\n", + " .add_horizon(zd=70, color=\"red\")\n", + " .add_hovertext()\n", + ")\n", + "viewable = builder.build()\n", + "bokeh.io.show(viewable)" + ] + }, + { + "cell_type": "markdown", + "id": "20", + "metadata": {}, + "source": [ + "## Sample with custom layout\n", + "\n", + "Simple adjustments to the layout can be made by passing a different class in the `bokeh.layout` module to the `layout` argument of `VisitMapBuilder.build`.\n", + "\n", + "More complex layout adjustments can be made by making subclass of `VisitMapBuilder` to repalce the `build` method with one that builds whatever layout is desired, for example:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "21", + "metadata": {}, + "outputs": [], + "source": [ + "class CustomVisitMapBuilder(schedview.plot.visit_skymaps.VisitMapBuilder):\n", + " def build(self) -> bokeh.models.UIElement:\n", + " # To have the plots that depend on them update when the\n", + " # alt, az, RA, Decl, or MJD sliders are adjusted, the\n", + " # `_connect_controls()` method needs to be called in\n", + " # the `build` method.\n", + " self._connect_controls()\n", + "\n", + " map_column_contents = [bokeh.models.Div(text=\"

Maps

\")]\n", + " control_column_contents = [bokeh.models.Div(text=\"

Controls

\")]\n", + "\n", + " # Some controls are associated with multiple plots.\n", + " # \"shown_control_names\" keepn track of which ones have been shown \n", + " # already on previous plots so we don't show shared controls\n", + " # multiple times.\n", + " shown_control_names = set()\n", + "\n", + " for sphere_map in self.spheremaps:\n", + " # force_update_time is an invisible div that\n", + " # needs to be included somewhere or bad things happen.\n", + " map_column_contents.append(sphere_map.force_update_time)\n", + "\n", + " # The .plot member of each SphereMap (or subclass) instance\n", + " # holds the instance of bokeh.plotting.figure (a subclass\n", + " # of bokeh.models.plots.Plot) that holds the figure itself.\n", + " map_column_contents.append(sphere_map.plot)\n", + "\n", + " # The bokeh models with the controls associated\n", + " # with a map are in a dictionar in the \"controls\"\n", + " # attribute.\n", + " for control_name in sphere_map.visible_control_names:\n", + " if control_name not in shown_control_names:\n", + " control_column_contents.append(sphere_map.controls[control_name])\n", + " shown_control_names.add(control_name)\n", + "\n", + " figure = bokeh.layouts.column(\n", + " [\n", + " bokeh.models.Div(text=\"

My custom figure output

\"),\n", + " bokeh.layouts.row(\n", + " [bokeh.layouts.column(map_column_contents), bokeh.layouts.column(control_column_contents)]\n", + " ),\n", + " ]\n", + " )\n", + " return figure\n", + "\n", + "\n", + "builder = (\n", + " CustomVisitMapBuilder(\n", + " visits,\n", + " mjd=visits[\"observationStartMJD\"].max(),\n", + " map_classes=[uranography.planisphere.Planisphere, uranography.armillary.ArmillarySphere],\n", + " )\n", + " .add_visit_patches()\n", + " .add_graticules()\n", + " .hide_future_visits()\n", + " .add_horizon()\n", + " .add_horizon(zd=70, color=\"red\")\n", + ")\n", + "\n", + "viewable = builder.build()\n", + "bokeh.io.show(viewable)" + ] + }, + { + "cell_type": "markdown", + "id": "22", + "metadata": {}, + "source": [ + "## Example with extra points and controls\n", + "\n", + "We can extend the `VisitMapBuilder` with new methods for new types of data or controls with new callbacks either by using `uranography` calls on the underlying `SphereMap` instances in the `spheremaps` attribute, or using the full `bokeh` API with the `plot` member of these `SphereMap` instances.\n", + "\n", + "For example, this example uses the plain `bokeh` API to add methods to mark the North and South celestial poles, and a toggle button to turn on and off the lables:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "23", + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd\n", + "import bokeh\n", + "import bokeh.models\n", + "\n", + "class CustomVisitMapBuilder(schedview.plot.visit_skymaps.VisitMapBuilder):\n", + "\n", + " def add_poles(self):\n", + " label_tuples = [('North Celestial Pole', 0.0, 90.0), ('South Celestial Pole', 0.0, -90.0)]\n", + " pole_data_source = bokeh.models.ColumnDataSource(pd.DataFrame(label_tuples, columns=['name', 'ra', 'decl']))\n", + "\n", + " # Save the pole bokeh objects to make it\n", + " # easier to connect callbacks later.\n", + " self.pole_labels = []\n", + "\n", + " # Add new elements to each plot in turn.\n", + " for spheremap in self.spheremaps:\n", + " # Add markers\n", + " spheremap.plot.scatter(\n", + " x=spheremap.proj_transform(\"x\", pole_data_source),\n", + " y=spheremap.proj_transform(\"y\", pole_data_source),\n", + " source=pole_data_source,\n", + " fill_color=\"lightblue\",\n", + " fill_alpha=0.0,\n", + " line_color=\"gray\",\n", + " marker=\"circle_cross\",\n", + " size=15\n", + " )\n", + "\n", + " # Add labels\n", + " labels = bokeh.models.LabelSet(\n", + " x=spheremap.proj_transform(\"x\", pole_data_source),\n", + " y=spheremap.proj_transform(\"y\", pole_data_source),\n", + " text=\"name\",\n", + " source=pole_data_source,\n", + " x_offset=6,\n", + " y_offset=6,\n", + " text_font_size=\"10pt\",\n", + " text_color=\"black\",\n", + " background_fill_color=\"white\",\n", + " background_fill_alpha=0.6,\n", + " border_line_color=None,\n", + " )\n", + " spheremap.plot.add_layout(labels)\n", + "\n", + " # Save the new labels to the list\n", + " # we will update with the callback.\n", + " self.pole_labels.append(labels)\n", + " \n", + "\n", + " # Make sure locations of the markers and labels\n", + " # get updated when the user changes the ra/dec or alt/az\n", + " # sliders.\n", + " spheremap.connect_controls(pole_data_source)\n", + "\n", + " return self\n", + "\n", + " def add_pole_label_toggle(self):\n", + " try:\n", + " pole_labels = self.pole_labels\n", + " except AttributeError:\n", + " raise UserError(\"add_poles must be called before add_pole_label_toggle\")\n", + "\n", + " # create a new control\n", + " toggle = bokeh.models.Toggle(label=\"Hide pole labels\", button_type=\"primary\", active=True)\n", + "\n", + " # Update the labels when the toggle state is changed.\n", + " js_code = \"\"\"\n", + " for (const label of labels) {\n", + " label.visible = toggle.active;\n", + " }\n", + " toggle.label = toggle.active ? \"Hide pole labels\" : \"Show pole labels\";\n", + " \"\"\"\n", + " for spheremap in self.spheremaps:\n", + " callback = bokeh.models.CustomJS(args=dict(labels=self.pole_labels, toggle=toggle), code=js_code)\n", + " toggle.js_on_change(\"active\", callback)\n", + "\n", + " # Add the new toggle to the list of things included\n", + " # shown in the layout with the reference map.\n", + " self.ref_map.controls['pole_label_toggle'] = toggle\n", + " \n", + " return self\n", + "\n", + "builder = (\n", + " CustomVisitMapBuilder(\n", + " visits,\n", + " mjd=visits[\"observationStartMJD\"].max(),\n", + " map_classes=[uranography.armillary.ArmillarySphere, uranography.planisphere.Planisphere],\n", + " figure_kwargs={\"match_aspect\": True},\n", + " )\n", + " .add_visit_patches()\n", + " .add_graticules()\n", + " .hide_future_visits()\n", + " .add_poles()\n", + " .add_pole_label_toggle()\n", + ")\n", + "viewable = builder.build()\n", + "bokeh.io.show(viewable)" + ] + }, + { + "cell_type": "markdown", + "id": "24", + "metadata": {}, + "source": [ + "## Example with customized colors" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "25", + "metadata": {}, + "outputs": [], + "source": [ + "applet_mode = False\n", + "\n", + "Planisphere = uranography.planisphere.Planisphere\n", + "ArmillarySphere = uranography.armillary.ArmillarySphere\n", + "VisitMapBuilder = schedview.plot.visit_skymaps.VisitMapBuilder\n", + "DARK_BAND_COLORS = schedview.plot.colors.LIGHT_PLOT_BAND_COLORS\n", + "\n", + "control_kwargs = {\"width\": None, \"styles\": {\"color\": \"#E5E5E5\"}}\n", + "# Planisphere and ArmillarySphere inherit default_slider_kwargs\n", + "# and default_select_kwargs from SphereMap, so setting these\n", + "# class attributes on SphereMap take care of both.\n", + "uranography.spheremap.SphereMap.default_slider_kwargs = control_kwargs\n", + "uranography.spheremap.SphereMap.default_select_kwargs = control_kwargs\n", + "\n", + "nside = 64\n", + "footprint_depth_by_band, footprint_regions = get_current_footprint(nside)\n", + "\n", + "maps = [ArmillarySphere, Planisphere] if not applet_mode else [Planisphere]\n", + "figure_specs = {\"match_aspect\": True, \"border_fill_color\": \"#262626\", \"background_fill_color\": \"#262626\"}\n", + "if applet_mode:\n", + " figure_specs.update({\"width\": 340, \"height\": 220})\n", + "\n", + "tooltips = \"\"\"\n", + "
\n", + "
Observation ID: @observationId
\n", + "
Start Timestamp: @start_timestamp{%F %T} UTC
\n", + "
Band: @band
\n", + "
RA, Dec: @fieldRA{0.000}, @fieldDec{0.000}
\n", + "
Observation Reason: @observation_reason
\n", + "
Science Program: @science_program
\n", + "
Para Angle: @paraAngle\\u00b0
\n", + "
azimulth, Altitude: @azimuth\\u00b0, @altitude\\u00b0
\n", + "
\n", + " \"\"\"\n", + "\n", + "builder = (\n", + " VisitMapBuilder(\n", + " visits,\n", + " mjd=visits[\"observationStartMJD\"].max(),\n", + " map_classes=maps,\n", + " visit_fill_colors=DARK_BAND_COLORS,\n", + " figure_kwargs=figure_specs,\n", + " )\n", + " .add_footprint_outlines(footprint_regions, line_width=5)\n", + " .add_visit_patches()\n", + " .hide_horizon_sliders()\n", + " .make_up_north()\n", + " .show_up_selector()\n", + " .add_eq_sliders()\n", + " .add_graticules()\n", + " .add_ecliptic()\n", + " .add_galactic_plane()\n", + " .add_datetime_slider()\n", + " .hide_mjd_slider()\n", + " .hide_future_and_other_night_visits()\n", + " .highlight_recent_visits()\n", + " .add_body(\"sun\", size=15, color=\"yellow\", alpha=1.0)\n", + " .add_body(\"moon\", size=15, color=\"orange\", alpha=0.8)\n", + " .add_horizon(color=\"#E5E5E5\", line_width=5)\n", + " .add_horizon(zd=70, color=\"red\", line_width=5)\n", + " .add_hovertext(visit_tooltips=tooltips)\n", + ")\n", + "\n", + "viewable = builder.build()\n", + "bokeh.io.show(viewable)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "26", + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "ehn_devel", + "language": "python", + "name": "ehn_devel" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.12.11" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/pyproject.toml b/pyproject.toml index 14fc1381..f58241c2 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -41,7 +41,7 @@ dependencies = [ "tabulate", "rubin-scheduler >= 3.22.0", "rubin-sim >= 2.6.1", - "uranography >= 1.4.1 ", + "uranography >= 1.5.0 ", ] [project.optional-dependencies] diff --git a/requirements.txt b/requirements.txt index 8b029d36..06222ca2 100644 --- a/requirements.txt +++ b/requirements.txt @@ -16,7 +16,7 @@ pytz scikit-learn skyproj tabulate -uranography>=1.4.1 +uranography>=1.5.0 rubin-scheduler>=3.22.0 rubin-sim>=2.6.1 pip diff --git a/schedview/compute/camera.py b/schedview/compute/camera.py index bc4cce9f..26c2f88c 100644 --- a/schedview/compute/camera.py +++ b/schedview/compute/camera.py @@ -2,6 +2,7 @@ import numpy as np import pandas as pd from astropy.coordinates import SkyCoord +from numpy.typing import NDArray class LsstCameraFootprintPerimeter(object): @@ -33,7 +34,7 @@ def __init__(self): self.vertices["angle"] = np.degrees(np.arctan2(self.vertices.y, self.vertices.x)) self.vertices["r"] = np.hypot(self.vertices.y, self.vertices.x) - def single_eq_vertices(self, ra, decl, rotation=0): + def single_eq_vertices(self, ra: float, decl: float, rotation: float = 0.0): """Compute vertices for a single pair of equatorial coordinates Parameters @@ -66,7 +67,12 @@ def single_eq_vertices(self, ra, decl, rotation=0): decl = eq_vertices.dec.deg return ra, decl - def __call__(self, ra, decl, rotation=0): + def __call__( + self, + ra: float | NDArray[np.floating], + decl: float | NDArray[np.floating], + rotation: float | NDArray[np.floating] = 0.0, + ): """Compute vertices for a single pair of equatorial coordinates Parameters @@ -88,13 +94,22 @@ def __call__(self, ra, decl, rotation=0): surrounding the camera footprints (degrees). """ if np.isscalar(ra): + assert isinstance(ra, float) + assert isinstance(decl, float) + assert isinstance(rotation, float) return self.single_eq_vertices(ra, decl, rotation) - if np.isscalar(rotation): - rotation = np.full_like(ra, rotation) + assert isinstance(ra, np.ndarray) or isinstance(ra, pd.Series) + assert isinstance(decl, np.ndarray) or isinstance(decl, pd.Series) + assert decl.shape == ra.shape + + rotation_array: NDArray[np.floating] | pd.Series = ( + np.full_like(ra, rotation) if np.isscalar(rotation) else rotation + ) + assert rotation_array.shape == ra.shape vertex_ras, vertex_decls = [], [] - for this_ra, this_decl, this_rotation in zip(ra, decl, rotation): + for this_ra, this_decl, this_rotation in zip(ra, decl, rotation_array): this_vertex_ra, this_vertex_decl = self.single_eq_vertices(this_ra, this_decl, this_rotation) vertex_ras.append(this_vertex_ra) vertex_decls.append(this_vertex_decl) diff --git a/schedview/plot/visit_skymaps.py b/schedview/plot/visit_skymaps.py new file mode 100644 index 00000000..8169d173 --- /dev/null +++ b/schedview/plot/visit_skymaps.py @@ -0,0 +1,1122 @@ +"""Utility for building visualizations of Rubin Observatory visits. + +Utility for creating interactive sky‑map visualisations of Rubin Observatory +observing visits using the **uranography** package. The module +defines builder class `VisitMapBuilder` following the "fluent builder" +pattern: a caller instantiates a `VisitMapBuilder` object, chains +configuration methods, and ultimately produces a `bokeh.models.UIElement` +that can be displayed in a report, dashboard, or other interface. +""" + +import warnings +from types import MethodType +from typing import Any, Callable, Dict, List, Optional, Self, SupportsFloat, Tuple + +import bokeh +import bokeh.layouts +import bokeh.models +import bokeh.palettes +import bokeh.plotting +import bokeh.transform +import healpy as hp +import numpy as np +import pandas as pd +from astropy.coordinates import SkyCoord, get_body +from astropy.time import Time +from uranography.api import ( + ArmillarySphere, + Planisphere, + SphereMap, + split_healpix_by_resolution, +) + +from schedview.collect import load_bright_stars +from schedview.compute.camera import LsstCameraFootprintPerimeter +from schedview.compute.footprint import find_healpix_area_polygons +from schedview.plot import PLOT_BAND_COLORS + +DEFAULT_VISIT_TOOLTIPS = ( + "@observationId: @start_timestamp{%F %T} UTC (mjd=@observationStartMJD{00000.0000}, " + + "LST=@observationStartLST\u00b0): " + + "@observation_reason (@science_program), " + + "in @band at \u03b1,\u03b4=@fieldRA\u00b0,@fieldDec\u00b0; " + + "q=@paraAngle\u00b0; a,A=@azimuth\u00b0,@altitude\u00b0" +) + +DEFAULT_VISIT_PATCHES_KWARGS = {"name": "visit_patches", "line_alpha": 0.0} + +SPHEREMAP_FIGURE_KWARGS = {"match_aspect": True} + + +class VisitMapBuilder: + """Builder for interactive visit sky‑maps. + + The ``VisitMapBuilder`` class follows the fluent builder pattern: + a ``VisitMapBuilder`` instance is created from a ``pandas.DataFrame`` of + Rubit Observatory visits. Callers then call series of chained + configuration methods to customise the visualisation. + The ``build`` method can then be used to return a `bokeh.models.UIElement` + instance that can be embedded in report, notebook, or dashboard. + + Parameters + ---------- + visits : `pandas.DataFrame` + Table of visits. At minimum the DataFrame must contain the columns + listed below. Columns can be renamed by overriding the + corresponding class attributes (e.g. ``VisitMapBuilder.ra_column``). + + Required columns + ^^^^^^^^^^^^^^^^ + ``fieldRA`` : ``float`` + Right‑ascension of the field pointing (degrees). + ``fieldDec`` : ``float`` + Declination of the field pointing (degrees). + ``observationStartMJD`` : ``float`` + Start time of the exposure as a Modified Julian Date. + ``band`` : ``str`` + Photometric band (e.g. ``u``, ``g``, ``r``, ``i``, ``z`` or ``y``). + ``rotSkyPos`` : ``float`` + Camera rotation angle (degrees). + + mjd : float, optional + Reference MJD for the underlying ``uranography`` maps. If omitted the + maximum ``observationStartMJD`` in *visits* is used; if *visits* is + empty the current time is used. + + map_classes : list, optional + List of ``uranography`` map classes to instantiate. By default an + :class:`~uranography.api.ArmillarySphere` and a + :class:`~uranography.api.Planisphere` are created. + + camera_perimeter : callable, optional + Function that returns camera‑footprint edge coordinates given arrays of + ``fieldRA``, ``fieldDec`` and ``rotSkyPos``. The default is + :class:`~schedview.compute.camera.LsstCameraFootprintPerimeter`. + + visit_fill_colors : Dict[str, str], optional + Colors for each of the bands. The default is None, which causes plots + to default to using ``schedview.plot.PLOT_BAND_COLORS``. + + mjd_slider_kwargs : Dict[str, Any], optional + Keyword arguments passed to `bokeh.models.Slider` for the MJD Slider. + + figure_kwargs : Dict[str, Asy], optional + Keword arguments with which to instantiate `bokeh.plotting.figure` + instances for each spheremap. The default is None, which causes plots + to default to using ``SPHEREMAP_FIGURE_KWARGS``. + + Notes + ----- + * The builder mutates its internal ``spheremaps`` list in‑place; each + ``add_`` method returns ``self`` to enable fluent chaining. + * Configuration methods such as `add_graticules` and `add_horizon` + are thin wrappers around the corresponding ``uranography`` API calls. + * Figure elements are rendered in the order in which the method calls + that add them are called, so that elements added later in the call + chain appear on top of elements that appear earlier. + + Examples + -------- + + >>> from uranography.api import ArmillarySphere, Planisphere + >>> import bokeh.io + >>> # Assume ``visits`` is a pandas.DataFrame with the required columns + >>> builder = ( + ... VisitMapBuilder( + ... visits, + ... mjd=visits['observationStartMJD'].max(), + ... map_classes=[ArmillarySphere, Planisphere] + ... ) + ... .add_visit_patches() + ... .add_graticules() + ... .add_ecliptic() + ... .add_galactic_plane() + ... .add_datetime_slider() + ... .add_eq_sliders() + ... .hide_future_visits() + ... .highlight_recent_visits() + ... .add_body('sun', size=15, color='yellow', alpha=1.0) + ... .add_body('moon', size=15, color='orange', alpha=0.8) + ... .add_horizon() + ... .add_horizon(zd=70, color='red') + ... ) + ... + >>> # Build the Bokeh layout + >>> viewable = builder.build() + >>> # Write the layout to a standalone HTML file + >>> bokeh.io.save(viewable, filename="visit_skymap.html") + """ + + mjd_column: str = "observationStartMJD" + ra_column: str = "fieldRA" + decl_column: str = "fieldDec" + rot_column: str = "rotSkyPos" + band_column: str = "band" + recent_max: float = 1.0 + recent_fade_scale: float = 2.0 / (24 * 60) + visit_columns: List[str] = [ + "observationId", + "start_timestamp", + "observationStartMJD", + "observationStartLST", + "band", + "fieldRA", + "fieldDec", + "rotSkyPos", + "paraAngle", + "azimuth", + "altitude", + "observation_reason", + "science_program", + ] + + def __init__( + self, + visits: pd.DataFrame | None = None, + mjd: Optional[float] = None, + map_classes: List[SphereMap] = [ArmillarySphere, Planisphere], + camera_perimeter: Optional[ + Callable[[np.ndarray, np.ndarray, np.ndarray], Tuple[np.ndarray, np.ndarray]] + ] = None, + visit_fill_colors: Optional[Dict[str, str]] = None, + mjd_slider_kwargs: Optional[Dict[str, Any]] = None, + figure_kwargs: Optional[Dict[str, Any]] = None, + ) -> None: + if mjd_slider_kwargs is None: + mjd_slider_kwargs = {} + + self.figure_kwargs: Dict[str, Any] = ( + SPHEREMAP_FIGURE_KWARGS if figure_kwargs is None else figure_kwargs + ) + self.visit_fill_colors: Dict[str, str] = ( + PLOT_BAND_COLORS if visit_fill_colors is None else visit_fill_colors + ) + self.visits = visits + + # If the mjd is not set by the caller, guess it from the visits + # if there are any, and if there are not, assume "now". + if mjd is None: + if self.visits is not None and len(self.visits) > 0: + self.mjd = self.visits[self.mjd_column].max() if mjd is None else mjd + else: + mjd_value = Time.now().mjd + # Tell the type checker that no, really, I know + # Time.now().mjd returns a scalar. + assert isinstance(mjd_value, SupportsFloat) + self.mjd = float(mjd_value) + else: + self.mjd = mjd + + self.instantiate_spheremaps(map_classes) + + self.camera_perimeter = ( + camera_perimeter if camera_perimeter is not None else LsstCameraFootprintPerimeter() + ) + + self.visits_ds = {} + self.visit_patches_added = False + self._add_mjd_slider(**mjd_slider_kwargs) + self.body_ds: Dict[str, bokeh.models.ColumnDataSource] = {} + self.horizon_ds: Dict[float, bokeh.models.ColumnDataSource] = {} + self.star_ds: Optional[bokeh.models.ColumnDataSource] = None + self.healpix_high_ds: Optional[bokeh.models.ColumnDataSource] = None + self.healpix_low_ds: Optional[bokeh.models.ColumnDataSource] = None + + def instantiate_spheremaps(self, map_classes: List[SphereMap]) -> None: + """Instantiate spheremap objects for each class in map_classes. + + This method creates spheremap instances for each class provided in the + `map_classes` list and stores them in `self.spheremaps`. The first + spheremap instance is also stored in `self.ref_map` for reference. + + Parameters + ---------- + map_classes : `list` [`uranography.api.SphereMap`] + A list of spheremap class constructors (e.g., ``ArmillarySphere``, + ``Planisphere``) to instantiate. Each class should accept a ``mjd`` + keyword argument. + + Notes + ----- + * This method is separated from `__init__` to allow subclasses to + override it and customize the spheremap instantiation process. + * The `self.spheremaps` list is populated with instances of the + provided classes, and `self.ref_map` is set to the first instance. + """ + # Separated into its own method so that it can be easily modified in + # subclasses. + # For example, a subclass could generate instances of bokeh figure to + # pass to each spheremap with customized parameters. + self.spheremaps = [ + mc(mjd=self.mjd, plot=bokeh.plotting.figure(**self.figure_kwargs)) for mc in map_classes + ] + self.ref_map = self.spheremaps[0] + + def add_visit_patches(self, visits: pd.DataFrame | None = None, **kwargs: Any) -> Self: + """Add visit patches to the map. + + Parameters + ---------- + visits : `pandas.DataFrame` or `None` + Table of visits. At minimum the DataFrame must contain the columns + listed below. Columns can be renamed by overriding the + corresponding class attributes + (e.g. ``VisitMapBuilder.ra_column``). + + Required columns + ^^^^^^^^^^^^^^^^ + ``fieldRA`` : ``float`` + Right‑ascension of the field pointing (degrees). + ``fieldDec`` : ``float`` + Declination of the field pointing (degrees). + ``observationStartMJD`` : ``float`` + Start time of the exposure as a Modified Julian Date. + ``band`` : ``str`` + Photometric band (e.g. ``u``, ``g``, ``r``, ``i``, ``z`` + or ``y``). + ``rotSkyPos`` : ``float`` + Camera rotation angle (degrees). + **kwargs + Additional keyword arguments passed to the underlying + `bokeh.plotting.figure.patches` call. + + Returns + ------- + self: `VisitMapBuilder` + Returns self to support method chaining. + """ + + if visits is not None: + if self.visits is not None: + raise ValueError("visits already supplied.") + self.visits = visits + + if self.visits is None: + warnings.warn("No visits to show.") + return self + + present_visit_columns = [c for c in self.visit_columns if c in self.visits.columns] + for band in "ugrizy": + in_band_mask = self.visits[self.band_column] == band + band_visits = self.visits.loc[in_band_mask, present_visit_columns].copy() + + if len(band_visits) < 1: + continue + + assert pd.api.types.is_float_dtype(band_visits[self.ra_column]) + assert pd.api.types.is_float_dtype(band_visits[self.decl_column]) + assert pd.api.types.is_float_dtype(band_visits[self.rot_column]) + # np.asarray makes type checking happier than .values + visits_ra = np.asarray(band_visits[self.ra_column]) + visits_decl = np.asarray(band_visits[self.decl_column]) + visits_rot = np.asarray(band_visits[self.rot_column]) + + ras, decls = self.camera_perimeter(visits_ra, visits_decl, visits_rot) + band_visits = band_visits.assign( + ra=ras, + decl=decls, + mjd=band_visits[self.mjd_column].values, + ) + + patches_kwargs = {"fill_color": self.visit_fill_colors[band]} + patches_kwargs.update(DEFAULT_VISIT_PATCHES_KWARGS) + patches_kwargs.update(kwargs) + + self.visits_ds[band] = self.ref_map.add_patches( + band_visits, + patches_kwargs=patches_kwargs, + ) + + for spheremap in self.spheremaps[1:]: + spheremap.add_patches(data_source=self.visits_ds[band], patches_kwargs=patches_kwargs) + + self.visit_patches_added = True + + return self + + def _add_mjd_slider(self, *args, **kwargs) -> Self: + """Add a slider to control the Date (MJD) of the visualization. + + This method adds an MJD slider to the visualization that allows users + to control the time displayed in the sky map. The slider is linked + across all spheremaps in the builder. + + Parameters + ---------- + *args : Any, optional + Positional arguments passed to `bokeh.models.Slider` + for the MJD Slider. + **kwargs : Dict[str, Any], optional + Keyword arguments passed to `bokeh.models.Slider` + for the MJD Slider. + + Returns + ------- + self: `VisitMapBuilder` + Returns self to support method chaining. + """ + mjd_now = Time.now().mjd + + # Appease type checker + assert isinstance(mjd_now, SupportsFloat) + + slider_kwargs = { + "start": self.visits[self.mjd_column].min() if self.visits is not None else float(mjd_now) - 1, + "end": self.visits[self.mjd_column].max() if self.visits is not None else float(mjd_now) + 1, + } + slider_kwargs.update(kwargs) + + if "mjd" not in self.ref_map.sliders: + self.ref_map.add_mjd_slider(*args, **slider_kwargs) + + self.mjd_slider = self.ref_map.sliders["mjd"] + + for spheremap in self.spheremaps[1:]: + spheremap.controls["mjd"] = self.mjd_slider + spheremap.suppressed_controls.append("mjd") + + # Support method chaining + return self + + def hide_mjd_slider(self) -> Self: + """Hide any mjd sliders. + + Returns + ------- + self: `VisitMapBuilder` + Returns self to support method chaining. + """ + for spheremap in self.spheremaps: + if "mjd" in spheremap.sliders: + spheremap.sliders["mjd"].visible = False + + return self + + def add_datetime_slider(self, *args: Any, **kwargs: Any) -> Self: + """Add a datetime slider linked to the (maybe invisible) MJD slider. + + Parameters + ---------- + *args + Additional positional arguments passed to the underlying + uranography datetime slider implementation. + **kwargs + Additional keyword arguments passed to the underlying uranography + datetime slider implementation. + + Returns + ------- + self : `VisitMapBuilder` + Returns self to support method chaining. + + Notes + ----- + + The datetime slider controls the same MJD value as the underlying MJD + slider, but displays it in a more traditional date/time format. + """ + self.ref_map.add_datetime_slider(*args, **kwargs) + return self + + def _hide_visits_by_mjd(self, hide_js: str, show_alpha: float = 0.5, hide_alpha: float = 0.0) -> Self: + """Hide visits based on MJD. + + Parameters + ---------- + hide_js : `str` + Javascript that shows or hides visits based on mjd. + show_alpha: `float` + The alpha of shown visits. + hide_alpha: `float` + The alpha of hidden visits. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable method chaining. + """ + + if "mjd" not in self.ref_map.sliders: + # The slider must exist for this feature to work, but if we + # have not yet explicitly added it, make it invisible. + self.ref_map.add_mjd_slider(visible=False) + + past_future_transform = bokeh.models.CustomJSTransform( + args=dict( + mjd_slider=self.ref_map.sliders["mjd"], + show_value=show_alpha, + hide_value=hide_alpha, + ), + v_func=hide_js, + ) + + # Apply the transform to the visit patches + try: + for spheremap in self.spheremaps: + visit_renderers = spheremap.plot.select(name="visit_patches") + if visit_renderers: + for renderer in visit_renderers: + renderer.glyph.fill_alpha = bokeh.transform.transform("mjd", past_future_transform) + except Exception as e: + warnings.warn(f"Could not apply hide_future_visits transform: {e}") + return self + + return self + + def hide_future_visits(self) -> Self: + """Hide visits that occur after the MJD value set by the slider. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable method chaining. + + Notes + ----- + * This method causes the map to show all visits prior to the MJD + set by the slider, regardless of what night they were observed in. + In contrast, the `VisitMapBuilder.hide_future_and_other_night_visits`, + hides both visits after MJD and also visits on other nights. + * If the MJD slider has not yet been added to the reference map, + this method adds an invisible MJD slider before creating the + transform. + """ + if not self.visit_patches_added: + raise ValueError("add_visit_patches must be called before hide_future_visits") + + # Transforms for recent, past, future visits + past_future_js = """ + const result = new Array(xs.length) + for (let i = 0; i < xs.length; i++) { + if (mjd_slider.value >= xs[i]) { + result[i] = show_value + } else { + result[i] = hide_value + } + } + return result + """ + + self._hide_visits_by_mjd(past_future_js) + return self + + def hide_future_and_other_night_visits(self) -> Self: + """Hide visits that occur after the MJD value set by the slider. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable method chaining. + + Notes + ----- + This method hides both visits that occurred on nights other than + that of MJD (including previous nights), and also visits that occur + later than MJD (even if in the same night). + In contrast, the `VisitMapBuilder.hide_future_visits`, + hides only visits that occur after MJD, and shows those that happened + on previous nights. + """ + # Transforms for recent, past, future visits + # The "- 0.5" in the floor causes it to use the same night rollover + # as dayobs. + past_future_js = """ + const result = new Array(xs.length) + for (let i = 0; i < xs.length; i++) { + if ((mjd_slider.value >= xs[i]) + && (Math.floor(mjd_slider.value - 0.5) == Math.floor(xs[i] - 0.5))) { + result[i] = show_value + } else { + result[i] = hide_value + } + } + return result + """ + + self._hide_visits_by_mjd(past_future_js) + return self + + def highlight_recent_visits(self) -> Self: + """Highlight recent visits based on the MJD slider. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable fluent method chaining. + + Notes + ----- + While the `VisitMapBuilder.hide_future_visits` and + `VisitMapBuilder.hide_future_and_other_night_visits` hide and show + visits entirely, this method highlights visits that occurred just + before MJD by giving them an extra outline. + """ + recent_js = """ + const result = new Array(xs.length) + for (let i = 0; i < xs.length; i++) { + if (mjd_slider.value < xs[i]) { + result[i] = 0 + } else { + result[i] = Math.max(0, max_value * (1 - (mjd_slider.value - xs[i]) / scale)) + } + } + return result + """ + + if not self.visit_patches_added: + raise ValueError("add_visit_patches must be called before highlight_recent_visits") + + if "mjd" not in self.ref_map.sliders: + # The slider must exist for this feature to work, but if we + # have not yet explicitly added it, make it invisible. + self.ref_map.add_mjd_slider(visible=False) + + recent_transform = bokeh.models.CustomJSTransform( + args=dict( + mjd_slider=self.ref_map.sliders["mjd"], + max_value=self.recent_max, + scale=self.recent_fade_scale, + ), + v_func=recent_js, + ) + + # Apply the transform to the visit patches + try: + for spheremap in self.spheremaps: + visit_renderers = spheremap.plot.select(name="visit_patches") + if visit_renderers: + for renderer in visit_renderers: + renderer.glyph.line_alpha = bokeh.transform.transform("mjd", recent_transform) + renderer.glyph.line_color = "#ff00ff" + renderer.glyph.line_width = 2 + except Exception as e: + warnings.warn(f"Could not apply highlight_recent_visits transform: {e}") + return self + + return self + + def decorate(self) -> Self: + """Add default decorations (graticules, ecliptic, and GP) to all maps. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable fluent method chaining. + """ + for spheremap in self.spheremaps: + spheremap.decorate() + + return self + + def add_ecliptic(self, *args: Any, **kwargs: Any) -> Self: + """Add the ecliptic to all maps. + + Parameters + ---------- + *args + Additional positional arguments forwarded to the underlying + ``add_ecliptic`` method of each ``SphereMap``. + **kwargs + Additional keyword arguments forwarded to the underlying + ``add_ecliptic`` method of each ``SphereMap``. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable fluent method chaining. + """ + for spheremap in self.spheremaps: + spheremap.add_ecliptic(*args, **kwargs) + + return self + + def add_galactic_plane(self, *args: Any, **kwargs: Any) -> Self: + """Add the galactic plane to all maps. + + Parameters + ---------- + *args + Additional positional arguments forwarded to the underlying + ``add_galactic_plane`` method of each ``SphereMap``. + **kwargs + Additional keyword arguments forwarded to the underlying + ``add_galactic_plane`` method of each ``SphereMap``. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable fluent method chaining. + """ + for spheremap in self.spheremaps: + spheremap.add_galactic_plane(*args, **kwargs) + + return self + + def add_graticules(self, *args: Any, **kwargs: Any) -> Self: + """Add graticules (RA/Dec grid lines) to all maps. + + Parameters + ---------- + *args + Additional positional arguments forwarded to the underlying + ``add_graticules`` method of each ``SphereMap``. + **kwargs + Additional keyword arguments forwarded to the underlying + ``add_graticules`` method of each ``SphereMap``. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable fluent method chaining. + """ + for spheremap in self.spheremaps: + spheremap.add_graticules(*args, **kwargs) + + return self + + def add_body(self, body: str, size: float, color: str, alpha: float, time_step: float = 1.0) -> Self: + """Add a celestial‑body marker to the reference map. + + Parameters + ---------- + body : `str` + Name of the solar system body (e.g. ``'sun'`` or ``'moon'``). + size : `float` + Marker size in screen pixels. + color : `str` + Color of the marker (any valid Bokeh color string). + alpha : `float` + Opacity of the marker (0.0‑1.0). + time_step: `float` + The time between "updates" to the position with time, in days. + + Returns + ------- + self : `VisitMapBuilder` + The builder instance to allow method chaining. + """ + + hide_js = """ + const result = new Array(xs.length) + for (let i = 0; i < xs.length; i++) { + if ((mjd_slider.value >= min_mjd) && (mjd_slider.value < max_mjd)) { + result[i] = show_value + } else { + result[i] = hide_value + } + } + return result + """ + + if time_step > self.mjd_slider.end - self.mjd_slider.start: + mjds = [(self.mjd_slider.end + self.mjd_slider.start) / 2] + else: + start_mjd: float = self.mjd_slider.start - time_step + end_mjd: float = self.mjd_slider.end + time_step + first_mjd: float = start_mjd + time_step / 2 + last_mjd: float = end_mjd + time_step / 2 + mjds = np.arange(first_mjd, last_mjd, time_step) + + ap_times = Time(mjds, format="mjd", scale="utc") + body_coords_all_times = get_body(body, ap_times) + assert isinstance(body_coords_all_times, SkyCoord) + + for mjd, ap_time, body_coords in zip(mjds, ap_times, body_coords_all_times): + assert isinstance(body_coords, SkyCoord) + assert isinstance(ap_time, Time) + + body_name = body if len(mjds) == 1 else body + ap_time.strftime("%Y%m%d%H%M%S") + min_mjd = mjd - time_step / 2 + max_mjd = mjd + time_step / 2 + + hide_js_transform = bokeh.models.CustomJSTransform( + args=dict( + mjd_slider=self.ref_map.sliders["mjd"], + min_mjd=min_mjd, + max_mjd=max_mjd, + show_value=alpha, + hide_value=0, + ), + v_func=hide_js, + ) + + hide_transform = bokeh.transform.transform("decl", hide_js_transform) + + circle_kwargs = { + "color": color, + "alpha": hide_transform, + "name": body, + } + + self.body_ds[body_name] = self.ref_map.add_marker( + ra=body_coords.ra.deg, + decl=body_coords.dec.deg, + name=body_name, + glyph_size=size, + circle_kwargs=circle_kwargs, + ) + + for spheremap in self.spheremaps[1:]: + spheremap.add_marker( + data_source=self.body_ds[body_name], + name=body_name, + glyph_size=size, + circle_kwargs=circle_kwargs, + ) + + return self + + def add_horizon(self, zd: float = 90, **line_kwargs: Any) -> Self: + """Add a horizon line to all maps. + + Parameters + ---------- + zd : `float`, optional + The zenith distance of the horizon line in degrees. Default is 90 + degrees (the horizon). + **line_kwargs + Additional keyword arguments passed to the underlying + ``add_horizon`` method of each ``SphereMap``. + + Returns + ------- + self : `VisitMapBuilder` + Returns ``self`` to enable fluent method chaining. + """ + data_source = self.ref_map.add_horizon(zd=zd, line_kwargs=line_kwargs) + for spheremap in self.spheremaps[1:]: + spheremap.add_horizon(data_source=data_source, line_kwargs=line_kwargs) + spheremap.connect_controls(data_source) + + self.horizon_ds[zd] = data_source + + return self + + def add_stars(self, star_data: Optional[pd.DataFrame] = None) -> Self: + """Add star markers to all maps. + + Parameters + ---------- + star_data : `pandas.DataFrame`, optional + DataFrame containing star data with columns "name", "ra", "decl", + and "Vmag". + If None, bright stars are loaded using `load_bright_stars()`. + Default is None. + + Returns + ------- + self : `VisitMapBuilder` + Returns self to enable method chaining. + """ + if star_data is None: + star_data = load_bright_stars().loc[:, ["name", "ra", "decl", "Vmag"]] + assert isinstance(star_data, pd.DataFrame) + + star_data["glyph_size"] = 15 - (15.0 / 3.5) * star_data["Vmag"] + star_data.query("glyph_size>0", inplace=True) + self.star_ds = self.ref_map.add_stars( + star_data, mag_limit_slider=False, star_kwargs={"color": "yellow"} + ) + + for spheremap in self.spheremaps[1:]: + spheremap.add_stars( + star_data, + data_source=self.star_ds, + mag_limit_slider=True, + star_kwargs={"color": "yellow"}, + ) + return self + + def add_hovertext(self, visit_tooltips: Optional[str] = None) -> Self: + """Add hover tooltips to visit patches. + + Parameters + ---------- + visit_tooltips : `str`, optional + The tooltip format string. If None, the default tooltip format + defined in DEFAULT_VISIT_TOOLTIPS is used. Default is None. + + Returns + ------- + self : `VisitMapBuilder` + Returns self to enable method chaining. + """ + if visit_tooltips is None: + visit_tooltips = DEFAULT_VISIT_TOOLTIPS + + for spheremap in self.spheremaps: + hover_tool = bokeh.models.HoverTool() + hover_tool.renderers = list(spheremap.plot.select({"name": "visit_patches"})) + hover_tool.tooltips = visit_tooltips + hover_tool.formatters = {"@start_timestamp": "datetime"} + spheremap.plot.add_tools(hover_tool) + return self + + def add_footprint(self, footprint: np.ndarray, nside_low: int = 8) -> Self: + """Add the LSST survey footprint to the sky maps. + + Parameters + ---------- + footprint : `numpy.ndarray` + A healpix map of the footprint. + nside_low : `int`, optional + The nside value for the low resolution component. Default is 8. + + Returns + ------- + self : `VisitMapBuilder` + Returns self to enable method chaining. + + Notes + ----- + Full healpix maps can make the full plots very large, resulting in + long transfer times and large files. + """ + cmap = bokeh.transform.linear_cmap("value", "Greys256", int(np.ceil(np.nanmax(footprint) * 2)), 0) + + nside_high = hp.npix2nside(footprint.shape[0]) + footprint_high, footprint_low = split_healpix_by_resolution(footprint, nside_low, nside_high) + + self.healpix_high_ds, cmap, glyph = self.ref_map.add_healpix( + footprint_high, nside=nside_high, cmap=cmap, name="footprint_high" + ) + self.healpix_low_ds, cmap, glyph = self.ref_map.add_healpix( + footprint_low, nside=nside_low, cmap=cmap, name="footprint_low" + ) + + for spheremap in self.spheremaps[1:]: + spheremap.add_healpix(self.healpix_high_ds, nside=nside_high, cmap=cmap, name="footprint_high") + spheremap.add_healpix(self.healpix_low_ds, nside=nside_low, cmap=cmap, name="footprint_low") + + for spheremap in self.spheremaps: + spheremap.connect_controls(self.healpix_high_ds) + spheremap.connect_controls(self.healpix_low_ds) + + return self + + @staticmethod + def _compute_footprint_outlines(footprint: np.ndarray) -> pd.DataFrame: + footprint_regions = footprint.copy() + footprint_regions[np.isin(footprint_regions, ["bulgy", "lowdust"])] = "WFD" + footprint_regions[ + np.isin(footprint_regions, ["LMC_SMC", "dusty_plane", "euclid_overlap", "nes", "scp", "virgo"]) + ] = "other" + + # Get rid of tiny little loops + footprint_outline = find_healpix_area_polygons(footprint_regions) + tiny_loops = footprint_outline.groupby(["region", "loop"]).count().query("RA<10").index + footprint_outline = footprint_outline.drop(tiny_loops) + return footprint_outline + + def add_footprint_outlines( + self, + footprint: np.ndarray, + colormap: Optional[Dict[str, str]] = None, + filled: bool = False, + **line_kwargs: Any, + ) -> Self: + """Add footprint outlines to the sky maps. + + Parameters + ---------- + footprint : `numpy.ndarray` + A HEALPix map of the survey footprint, where each pixel contains + a region identifier. + colormap : `dict` [`str`, `str`], optional + Dictionary mapping region names to colors. If None, a default + color palette is automatically generated based on the number of + unique regions in the footprint. Default is None. + filled : `bool`, optional + If True, fill the regions with color instead of drawing outlines. + If False, draw only the outline of each region. Default is False. + **line_kwargs + Additional keyword arguments passed to the underlying Bokeh + plotting methods (``plot.patch`` for filled regions or + ``plot.line`` for outlines). These can include line width, alpha, + etc. + + Returns + ------- + self : `VisitMapBuilder` + Returns self to enable method chaining. + + Notes + ----- + * When `filled=True`, the filled regions may not render correctly if + the polygon crosses a projection discontinuity. + * The footprint regions are defined by the values in the input + `footprint` array. Special regions like "bulgy", "lowdust", etc. are + grouped into "WFD" (Wide, Deep, and Fast), and other regions are + grouped into "other". + """ + footprint_polygons = self._compute_footprint_outlines(footprint) + + if filled: + warnings.warn( + 'The "filled" option does work correctly when the polygon crosses a projection discontinuity.' + ) + + if "region" not in footprint_polygons.index.names: + footprint_polygons = ( + footprint_polygons.assign(region="only").set_index("region", append=True).copy() + ) + + if "loop" not in footprint_polygons.index.names: + footprint_polygons = footprint_polygons.assign(loop=0).set_index("loop", append=True).copy() + + footprint_polygons = footprint_polygons.reorder_levels(["region", "loop"]).copy() + + outside = "" + if colormap is None: + regions = [ + r for r in footprint_polygons.index.get_level_values("region").unique() if r != outside + ] + if len(regions) == 1: + palette = ["black"] + elif len(regions) == 2: + palette = ["black", "darkgray"] + else: + # Try palettes from the "colorblind" section in the bokeh docs + try: + palette = bokeh.palettes.Colorblind[len(regions)] + except KeyError: + palette = bokeh.palettes.TolRainbow[len(regions)] + + colormap = {r: c for r, c in zip(regions, palette)} + + footprint_regions = {} + for region_index in set(footprint_polygons.index.values.tolist()): + assert isinstance(region_index, tuple) + region_name, loop_id = region_index + if region_name not in footprint_regions: + footprint_regions[region_name] = {} + footprint_regions[region_name][loop_id] = ( + footprint_polygons.loc[(region_name, loop_id), ["RA", "decl"]] + .apply(tuple, axis="columns") + .tolist() + ) + + for region_name in footprint_regions: + if region_name == outside: + continue + + for this_loop in footprint_regions[region_name].values(): + line_kwargs["name"] = "footprint_outline" + loop_ds = bokeh.models.ColumnDataSource({"coords": this_loop}) + for spheremap in self.spheremaps: + if filled: + spheremap.plot.patch( + spheremap.x_transform("coords"), + spheremap.y_transform("coords"), + color=colormap[region_name], + source=loop_ds, + **line_kwargs, + ) + else: + spheremap.plot.line( + spheremap.x_transform("coords"), + spheremap.y_transform("coords"), + color=colormap[region_name], + source=loop_ds, + **line_kwargs, + ) + spheremap.connect_controls(loop_ds) + + return self + + def add_eq_sliders(self) -> Self: + """Add sliders for setting the center using equatorial coordinates. + + Returns + ------- + self: `VisitMapBuilder` + Returns self to support method chaining. + """ + + for spheremap in self.spheremaps: + # Only call on maps for which it applies + maybe_add_eq_sliders_method = getattr(spheremap, "add_eq_sliders", None) + if isinstance(maybe_add_eq_sliders_method, MethodType): + spheremap.add_eq_sliders() + + return self + + def hide_horizon_sliders(self) -> Self: + """Set horizon sliders to not be visible. + + Returns + ------- + self: `VisitMapBuilder` + Returns self to support method chaining. + """ + for spheremap in self.spheremaps: + for coord in ("alt", "az"): + if coord in spheremap.sliders: + spheremap.sliders[coord].visible = False + + return self + + def make_up_north(self) -> Self: + """Set "up" in the maps to be north. + + Returns + ------- + self: `VisitMapBuilder` + Returns self to support method chaining. + """ + + for spheremap in self.spheremaps: + if "up" in spheremap.sliders: + spheremap.sliders["up"].value = "north is up" + + return self + + def show_up_selector(self) -> Self: + """Make the selector for the orientation visible. + + Returns + ------- + self: `VisitMapBuilder` + Returns self to support method chaining. + """ + + for spheremap in self.spheremaps: + if "up" in spheremap.sliders: + spheremap.sliders["up"].visible = True + + return self + + def hide_up_selector(self) -> Self: + """Make the selector for the orientation invisible. + + Returns + ------- + self: `VisitMapBuilder` + Returns self to support method chaining. + """ + + for spheremap in self.spheremaps: + if "up" in spheremap.sliders: + spheremap.sliders["up"].visible = False + + return self + + def _connect_controls(self): + # Must be called after all data sources and contros have been added. + + # Data sources that might change with silders on all plots, not just + # those that change with orientation. + dynamic_data_sources = list(self.visits_ds.values()) + list(self.body_ds.values()) + for spheremap in self.spheremaps: + for data_source in dynamic_data_sources: + spheremap.connect_controls(data_source) + + def build( + self, layout: Callable[[List[bokeh.models.UIElement]], bokeh.models.UIElement] = bokeh.layouts.row + ) -> bokeh.models.UIElement: + + self._connect_controls() + + map_figures = list(s.figure for s in self.spheremaps) + combined_figure = layout(map_figures) + return combined_figure diff --git a/schedview/plot/visitmap.py b/schedview/plot/visitmap.py index 472f708d..298c376e 100644 --- a/schedview/plot/visitmap.py +++ b/schedview/plot/visitmap.py @@ -347,6 +347,12 @@ def create_visit_skymaps( The arguments used to produce the figure using `plot_visit_skymaps`. """ + warn( + "create_visit_skymaps is deprecated and will be removed; use VisitMapBuilder instead.", + DeprecationWarning, + stacklevel=2, + ) + site = None if observatory is None else observatory.location night_events = schedview.compute.astro.night_events(night_date=night_date, site=site, timezone=timezone) start_time = Time(night_events.loc["sunset", "UTC"]) diff --git a/tests/test_plot_visit_skymaps.py b/tests/test_plot_visit_skymaps.py new file mode 100644 index 00000000..6e9daafc --- /dev/null +++ b/tests/test_plot_visit_skymaps.py @@ -0,0 +1,237 @@ +"""Unit tests for schedview.plot.visit_skymaps module.""" + +import pathlib +import tempfile + +import bokeh.io +import bokeh.models +import pandas as pd +from rubin_scheduler.scheduler.utils import get_current_footprint +from uranography.api import ArmillarySphere, Planisphere + +from schedview.plot.visit_skymaps import VisitMapBuilder + +TEST_VISITS = pd.DataFrame( + { + "fieldRA": [10.0, 20.0], + "fieldDec": [-5.0, 15.0], + "observationStartMJD": [59000.0, 59001.0], + "band": ["g", "r"], + "rotSkyPos": [0.0, 45.0], + "observationId": [1, 2], + "start_timestamp": pd.to_datetime(["2020-01-01", "2020-01-02"]), + "observationStartLST": [150.0, 160.0], + "paraAngle": [0.0, 0.0], + "azimuth": [180.0, 190.0], + "altitude": [45.0, 50.0], + "observation_reason": ["science", "science"], + "science_program": ["prog1", "prog2"], + } +) +FOOTPRINT_NSIDE = 16 + + +def _save_and_check_viewable_html( + viewable: bokeh.models.UIElement, filename: str = "visit_skymap.html" +) -> None: + """ + Save a Bokeh viewable to a temporary HTML file and assert that the file + was written and has non-zero size. + + Parameters + ---------- + viewable: bokeh.models.UIElement + The Bokeh object returned by ``VisitMapBuilder.build()``. + filename: str, optional + Name of the file inside the temporary directory (default + ``"visit_skymap.html"``). + """ + with tempfile.TemporaryDirectory() as temp_dir: + out_file = pathlib.Path(temp_dir) / filename + bokeh.io.save(viewable, filename=str(out_file)) + assert out_file.is_file() + assert out_file.stat().st_size > 0 + + +def test_basic_build(): + """Test that the most basic builder builds.""" + builder = VisitMapBuilder(TEST_VISITS) + builder.add_visit_patches() + viewable = builder.build() + + assert isinstance(viewable, bokeh.models.UIElement) + + visit_patches = viewable.select({"name": "visit_patches"}) + assert len(list(visit_patches)) > 0 + + _save_and_check_viewable_html(viewable) + + +def test_elaborate_build(): + """Test a builder with many chained options""" + visits = TEST_VISITS + footprint_regions = get_current_footprint(FOOTPRINT_NSIDE)[1] + builder = ( + VisitMapBuilder( + visits, mjd=visits["observationStartMJD"].max(), map_classes=[ArmillarySphere, Planisphere] + ) + .add_visit_patches() + .add_graticules() + .add_ecliptic() + .add_galactic_plane() + ._add_mjd_slider(start=visits["observationStartMJD"].min(), end=visits["observationStartMJD"].max()) + .add_datetime_slider() + .hide_future_visits() + .highlight_recent_visits() + .add_footprint_outlines(footprint_regions) + .add_body("sun", size=15, color="yellow", alpha=1.0) + .add_body("moon", size=15, color="orange", alpha=0.8) + .add_horizon() + .add_horizon(zd=70, color="red") + .add_eq_sliders() + .hide_horizon_sliders() + .make_up_north() + .show_up_selector() + ) + viewable = builder.build() + assert isinstance(viewable, bokeh.models.UIElement) + _save_and_check_viewable_html(viewable) + + +def test_add_graticules(): + """Test that add_graticules adds graticule renderers.""" + builder = VisitMapBuilder(TEST_VISITS) + builder.add_visit_patches() + builder.add_graticules() + viewable = builder.build() + + graticule_glyphs = viewable.select({"name": "graticule_glyph"}) + assert len(list(graticule_glyphs)) > 0 + + _save_and_check_viewable_html(viewable) + + +def test_add_ecliptic(): + """Test that add_ecliptic adds ecliptic renderers.""" + builder = VisitMapBuilder(TEST_VISITS) + builder.add_ecliptic() + viewable = builder.build() + + ecliptic_glyphs = viewable.select({"name": "ecliptic_glyph"}) + assert len(list(ecliptic_glyphs)) > 0 + + _save_and_check_viewable_html(viewable) + + +def test_add_body_sun(): + """Test that add_body adds a sun marker.""" + builder = VisitMapBuilder(TEST_VISITS) + builder.add_body("sun", size=10, color="yellow", alpha=1.0) + viewable = builder.build() + + sun_renderers = viewable.select({"name": "sun"}) + assert len(list(sun_renderers)) > 0 + + _save_and_check_viewable_html(viewable) + + +def test_add_horizon(): + """Test that add_horizon adds a horizon line.""" + builder = VisitMapBuilder(TEST_VISITS) + builder.add_horizon() + viewable = builder.build() + + horizon_glyphs = viewable.select({"name": "horizon_glyph"}) + assert len(list(horizon_glyphs)) > 0 + + _save_and_check_viewable_html(viewable) + + +def test_add_mjd_slider(): + """Test that add_mjd_slider initialises the mjd_slider attribute.""" + builder = VisitMapBuilder(TEST_VISITS) + builder._add_mjd_slider() + viewable = builder.build() + + assert isinstance(builder.mjd_slider, bokeh.models.Slider) + + _save_and_check_viewable_html(viewable) + + +def test_hide_future_visits(): + """Test that hide_future_visits applies a transform to visit patches.""" + builder = VisitMapBuilder(TEST_VISITS) + builder.add_visit_patches() + builder.hide_future_visits() + viewable = builder.build() + + visit_renderers = viewable.select({"name": "visit_patches"}) + assert len(list(visit_renderers)) > 0 + # TODO look for callback + + _save_and_check_viewable_html(viewable) + + +def test_highlight_recent_visits(): + """Test that highlight_recent_visits applies its transform.""" + builder = VisitMapBuilder(TEST_VISITS) + builder.add_visit_patches() + builder.highlight_recent_visits() + viewable = builder.build() + + # Check that a transform was applied to the visit patches + visit_renderers = viewable.select({"name": "visit_patches"}) + assert len(list(visit_renderers)) > 0 + # TODO look for callback + + _save_and_check_viewable_html(viewable) + + +def test_add_hovertext(): + """Test that add_hovertext attaches a HoverTool to visit patches.""" + builder = VisitMapBuilder(TEST_VISITS) + builder.add_hovertext() + viewable = builder.build() + + hover_tools = list(viewable.select({"type": bokeh.models.HoverTool})) + assert len(hover_tools) > 0 + + _save_and_check_viewable_html(viewable) + + +def test_add_footprint(): + footprint_depth_by_band = get_current_footprint(FOOTPRINT_NSIDE)[0] + bands_in_footprint = tuple(footprint_depth_by_band.dtype.fields.keys()) + footprint_depth = footprint_depth_by_band[bands_in_footprint[0]] + for band in bands_in_footprint[1:]: + footprint_depth = footprint_depth + footprint_depth_by_band[band] + + builder = VisitMapBuilder(TEST_VISITS) + builder.add_footprint(footprint_depth) + viewable = builder.build() + + for rend_name in ("footprint_high", "footprint_low"): + footprint_renderers = list(viewable.select({"name": rend_name})) + assert len(footprint_renderers) > 0 + + _save_and_check_viewable_html(viewable) + + +def test_add_footprint_outlines(): + footprint_regions = get_current_footprint(FOOTPRINT_NSIDE)[1] + + builder = VisitMapBuilder(TEST_VISITS) + builder.add_footprint_outlines(footprint_regions) + viewable = builder.build() + + # Bokeh does not support pattern matching in selection by name, + # so iterate over all renderers and check their names explicitly. + outline_renderers = [] + for renderer in viewable.select({"type": bokeh.models.GlyphRenderer}): + assert isinstance(renderer.name, str) + if renderer.name.startswith("footprint_outline"): + outline_renderers.append(renderer) + + assert len(outline_renderers) > 0 + + _save_and_check_viewable_html(viewable)