diff --git a/.github/workflows/wipac-cicd.yml b/.github/workflows/wipac-cicd.yml index 5ff32407..3db893e3 100644 --- a/.github/workflows/wipac-cicd.yml +++ b/.github/workflows/wipac-cicd.yml @@ -84,48 +84,6 @@ jobs: if: always() run: | more *.json | cat - - functionality-tests: - needs: [py-versions] - runs-on: ubuntu-latest - strategy: - fail-fast: false - matrix: - py3: ${{ fromJSON(needs.py-versions.outputs.matrix) }} - input-file: [ "millipede_wilks/run00136662.evt000035405932.neutrino_1" ] - steps: - - uses: actions/checkout@v3 - - uses: actions/setup-python@v4 - with: - python-version: ${{ matrix.py3 }} - - uses: actions/checkout@v3 - with: - repository: icecube/skymap_scanner - ref: main - path: skymap_scanner - sparse-checkout: tests/data - - name: pip install - run: | - pip install --upgrade pip wheel setuptools - pip install . - - name: Run conversion script - # working-directory: auxiliary-repo/path/to/directory # Replace 'path/to/directory' with the path to the directory containing your Python script - run: | - python examples/npz_to_json.py --input-file skymap_scanner/tests/data/results_npz/${{ matrix.input-file }}.npz --output-dir . - cat $(ls -Art | tail -n 1) # cat outfile - - # Wait for skymap_scanner PR#163 - # - name: Compare output - # run: | - # tree - # outfile=$(ls -Art | tail -n 1) - # diff $outfile skymap_scanner/tests/data/results_json/${{ matrix.input-file }}.json - - - - - - release: # only run on main/master/default diff --git a/examples/plot_skydriver_scan_result.py b/examples/plot_skydriver_scan_result.py index fd359d40..dc92a353 100644 --- a/examples/plot_skydriver_scan_result.py +++ b/examples/plot_skydriver_scan_result.py @@ -39,8 +39,8 @@ def main() -> None: serialzed = rc.request_seq("GET", f"/scan/{args.scan_id}/result")["skyscan_result"] result = SkyScanResult.deserialize(serialzed) - result.create_plot(dosave=True) - result.create_plot_zoomed(dosave=True, plot_bounding_box=True) + result.create_plot() + result.create_plot_zoomed(plot_bounding_box=True) if __name__ == "__main__": diff --git a/requirements-examples.txt b/requirements-examples.txt index 9c0c6128..77521e09 100644 --- a/requirements-examples.txt +++ b/requirements-examples.txt @@ -4,7 +4,7 @@ # # pip-compile --extra=examples --output-file=requirements-examples.txt # -astropy==5.3.1 +astropy==5.3.4 # via # healpy # icecube-skyreader (setup.py) @@ -12,31 +12,31 @@ cachetools==5.3.1 # via wipac-rest-tools certifi==2023.7.22 # via requests -cffi==1.15.1 +cffi==1.16.0 # via cryptography -charset-normalizer==3.2.0 +charset-normalizer==3.3.0 # via requests -contourpy==1.1.0 +contourpy==1.1.1 # via matplotlib -cryptography==41.0.3 +cryptography==41.0.4 # via pyjwt -cycler==0.11.0 +cycler==0.12.1 # via matplotlib -fonttools==4.42.0 +fonttools==4.43.1 # via matplotlib -healpy==1.16.3 +healpy==1.16.6 # via icecube-skyreader (setup.py) idna==3.4 # via requests -kiwisolver==1.4.4 +kiwisolver==1.4.5 # via matplotlib -matplotlib==3.7.2 +matplotlib==3.8.0 # via # healpy # icecube-skyreader (setup.py) meander==0.0.3 # via icecube-skyreader (setup.py) -numpy==1.25.2 +numpy==1.26.0 # via # astropy # contourpy @@ -46,13 +46,13 @@ numpy==1.25.2 # pandas # pyerfa # scipy -packaging==23.1 +packaging==23.2 # via # astropy # matplotlib -pandas==2.0.3 +pandas==2.1.1 # via icecube-skyreader (setup.py) -pillow==10.0.0 +pillow==10.0.1 # via matplotlib pycparser==2.21 # via cffi @@ -60,7 +60,7 @@ pyerfa==2.0.0.3 # via astropy pyjwt[crypto]==2.8.0 # via wipac-rest-tools -pyparsing==3.0.9 +pyparsing==3.1.1 # via matplotlib pypng==0.20220715.0 # via qrcode @@ -68,7 +68,7 @@ python-dateutil==2.8.2 # via # matplotlib # pandas -pytz==2023.3 +pytz==2023.3.post1 # via pandas pyyaml==6.0.1 # via astropy @@ -81,23 +81,25 @@ requests==2.31.0 # wipac-rest-tools requests-futures==1.0.1 # via wipac-rest-tools -scipy==1.11.1 +scipy==1.11.3 # via healpy six==1.16.0 # via python-dateutil -tornado==6.3.2 +tornado==6.3.3 # via wipac-rest-tools -typing-extensions==4.7.1 +typing-extensions==4.8.0 # via # qrcode # wipac-dev-tools tzdata==2023.3 # via pandas -urllib3==2.0.4 - # via requests -wipac-dev-tools==1.6.16 +urllib3==2.0.6 + # via + # requests + # wipac-rest-tools +wipac-dev-tools==1.7.0 # via # icecube-skyreader (setup.py) # wipac-rest-tools -wipac-rest-tools==1.5.1 +wipac-rest-tools==1.5.3 # via icecube-skyreader (setup.py) diff --git a/requirements-tests.txt b/requirements-tests.txt index 374604aa..0c0e244f 100644 --- a/requirements-tests.txt +++ b/requirements-tests.txt @@ -4,37 +4,37 @@ # # pip-compile --extra=tests --output-file=requirements-tests.txt # -astropy==5.3.1 +astropy==5.3.4 # via # healpy # icecube-skyreader (setup.py) certifi==2023.7.22 # via requests -charset-normalizer==3.2.0 +charset-normalizer==3.3.0 # via requests -contourpy==1.1.0 +contourpy==1.1.1 # via matplotlib -cycler==0.11.0 +cycler==0.12.1 # via matplotlib -exceptiongroup==1.1.2 +exceptiongroup==1.1.3 # via pytest -fonttools==4.42.0 +fonttools==4.43.1 # via matplotlib -healpy==1.16.3 +healpy==1.16.6 # via icecube-skyreader (setup.py) idna==3.4 # via requests iniconfig==2.0.0 # via pytest -kiwisolver==1.4.4 +kiwisolver==1.4.5 # via matplotlib -matplotlib==3.7.2 +matplotlib==3.8.0 # via # healpy # icecube-skyreader (setup.py) meander==0.0.3 # via icecube-skyreader (setup.py) -numpy==1.25.2 +numpy==1.26.0 # via # astropy # contourpy @@ -44,22 +44,22 @@ numpy==1.25.2 # pandas # pyerfa # scipy -packaging==23.1 +packaging==23.2 # via # astropy # matplotlib # pytest -pandas==2.0.3 +pandas==2.1.1 # via icecube-skyreader (setup.py) -pillow==10.0.0 +pillow==10.0.1 # via matplotlib -pluggy==1.2.0 +pluggy==1.3.0 # via pytest pyerfa==2.0.0.3 # via astropy -pyparsing==3.0.9 +pyparsing==3.1.1 # via matplotlib -pytest==7.4.0 +pytest==7.4.2 # via # icecube-skyreader (setup.py) # pytest-mock @@ -69,23 +69,23 @@ python-dateutil==2.8.2 # via # matplotlib # pandas -pytz==2023.3 +pytz==2023.3.post1 # via pandas pyyaml==6.0.1 # via astropy requests==2.31.0 # via wipac-dev-tools -scipy==1.11.1 +scipy==1.11.3 # via healpy six==1.16.0 # via python-dateutil tomli==2.0.1 # via pytest -typing-extensions==4.7.1 +typing-extensions==4.8.0 # via wipac-dev-tools tzdata==2023.3 # via pandas -urllib3==2.0.4 +urllib3==2.0.6 # via requests -wipac-dev-tools==1.6.16 +wipac-dev-tools==1.7.0 # via icecube-skyreader (setup.py) diff --git a/requirements.txt b/requirements.txt index 1df988c4..90a0beb6 100644 --- a/requirements.txt +++ b/requirements.txt @@ -4,33 +4,33 @@ # # pip-compile --output-file=requirements.txt # -astropy==5.3.1 +astropy==5.3.4 # via # healpy # icecube-skyreader (setup.py) certifi==2023.7.22 # via requests -charset-normalizer==3.2.0 +charset-normalizer==3.3.0 # via requests -contourpy==1.1.0 +contourpy==1.1.1 # via matplotlib -cycler==0.11.0 +cycler==0.12.1 # via matplotlib -fonttools==4.42.0 +fonttools==4.43.1 # via matplotlib -healpy==1.16.3 +healpy==1.16.6 # via icecube-skyreader (setup.py) idna==3.4 # via requests -kiwisolver==1.4.4 +kiwisolver==1.4.5 # via matplotlib -matplotlib==3.7.2 +matplotlib==3.8.0 # via # healpy # icecube-skyreader (setup.py) meander==0.0.3 # via icecube-skyreader (setup.py) -numpy==1.25.2 +numpy==1.26.0 # via # astropy # contourpy @@ -40,37 +40,37 @@ numpy==1.25.2 # pandas # pyerfa # scipy -packaging==23.1 +packaging==23.2 # via # astropy # matplotlib -pandas==2.0.3 +pandas==2.1.1 # via icecube-skyreader (setup.py) -pillow==10.0.0 +pillow==10.0.1 # via matplotlib pyerfa==2.0.0.3 # via astropy -pyparsing==3.0.9 +pyparsing==3.1.1 # via matplotlib python-dateutil==2.8.2 # via # matplotlib # pandas -pytz==2023.3 +pytz==2023.3.post1 # via pandas pyyaml==6.0.1 # via astropy requests==2.31.0 # via wipac-dev-tools -scipy==1.11.1 +scipy==1.11.3 # via healpy six==1.16.0 # via python-dateutil -typing-extensions==4.7.1 +typing-extensions==4.8.0 # via wipac-dev-tools tzdata==2023.3 # via pandas -urllib3==2.0.4 +urllib3==2.0.6 # via requests -wipac-dev-tools==1.6.16 +wipac-dev-tools==1.7.0 # via icecube-skyreader (setup.py) diff --git a/setup.cfg b/setup.cfg index ecaed838..83ca8da6 100644 --- a/setup.cfg +++ b/setup.cfg @@ -29,6 +29,7 @@ classifiers = Programming Language :: Python :: 3.9 Programming Language :: Python :: 3.10 Programming Language :: Python :: 3.11 + Programming Language :: Python :: 3.12 download_url = https://pypi.org/project/icecube-skyreader/ project_urls = Tracker = https://github.com/icecube/skyreader/issues @@ -53,7 +54,7 @@ install_requires = numpy pandas wipac-dev-tools -python_requires = >=3.8, <3.12 +python_requires = >=3.8, <3.13 packages = find: [options.extras_require] diff --git a/skyreader/plot/plotting_tools.py b/skyreader/plot/plotting_tools.py index 765859e5..377430ff 100644 --- a/skyreader/plot/plotting_tools.py +++ b/skyreader/plot/plotting_tools.py @@ -192,6 +192,7 @@ def set_xlim(self, *args, **kwargs): Axes.set_ylim(self, -np.pi / 2.0, np.pi / 2.0) def _get_core_transform(self, resolution): + # mypy error: "_get_core_transform" undefined in superclass [misc] return Affine2D().translate(-np.pi, 0.) + super(AstroMollweideAxes, self)._get_core_transform(resolution) class RaFormatter(Formatter): @@ -207,6 +208,7 @@ def __call__(self, x, pos=None): def set_longitude_grid(self, degrees): # Copied from matplotlib.geo.GeoAxes.set_longitude_grid and modified number = (360 // degrees) + 1 + # mypy error: Argument 1 to "FixedLocator" has incompatible type "ndarray[Any, dtype[floating[Any]]]"; expected "Sequence[float]" self.xaxis.set_major_locator( FixedLocator( np.linspace(0, 2*np.pi, number, True)[1:-1])) @@ -215,6 +217,7 @@ def set_longitude_grid(self, degrees): def _set_lim_and_transforms(self): # Copied from matplotlib.geo.GeoAxes._set_lim_and_transforms and modified + # mypy error: Argument 1 to "FixedLocator" has incompatible type "ndarray[Any, dtype[floating[Any]]]"; expected "Sequence[float]" super(AstroMollweideAxes, self)._set_lim_and_transforms() # This is the transform for latitude ticks. @@ -223,6 +226,8 @@ def _set_lim_and_transforms(self): self._yaxis_transform = \ yaxis_stretch + \ self.transData + # mypy error: "AstroMollweideAxes" has no attribute "transProjection" [attr-defined] + # mypy error: "AstroMollweideAxes" has no attribute "transAffine" [attr-defined] yaxis_text_base = \ yaxis_stretch + \ self.transProjection + \ diff --git a/skyreader/result.py b/skyreader/result.py index e5e0b98d..577df05f 100644 --- a/skyreader/result.py +++ b/skyreader/result.py @@ -4,8 +4,6 @@ # pylint: skip-file # flake8: noqa - -import io import itertools as it import json import logging @@ -584,47 +582,52 @@ def best_dir(self): Plotting routines """ - def create_plot(self, - dosave=False, - dozoom=False, - log_func=None, - upload_func=None, - final_channels=None): - - if log_func is None: - def log_func(x): - print(x) - - if upload_func is None: - def upload_func(file_buffer, name, title): - pass - - if final_channels is None: - final_channels=["#test_messaging"] - y_inches = 3.85 - x_inches = 6 - dpi = 150 if not dozoom else 1200 - xsize = x_inches*dpi - ysize = xsize//2 - - lonra=[-10.,10.] - latra=[-10.,10.] + plot_y_size_in = 3.85 + plot_x_size_in = 6 + plot_dpi_standard = 150 + plot_dpi_zoomed = 1200 + plot_colormap = matplotlib.colormaps['plasma_r'] + def check_result(self): + """Check in legacy plotting code. + """ for k in self.result: if "nside-" not in k: raise RuntimeError("\"nside\" not in result file..") + + @staticmethod + # Calculates are using Gauss-Green theorem / shoelace formula + # TODO: vectorize using numpy. + # Note: in some cases the argument is not a np.ndarray so one has to convert the data series beforehand. + def calculate_area(vs) -> float: + a = 0 + x0, y0 = vs[0] + for [x1,y1] in vs[1:]: + dx = x1-x0 + dy = y1-y0 + a += 0.5*(y0*dx - x0*dy) + x0 = x1 + y0 = y1 + return a + + def create_plot(self, dozoom = False): + + dpi = self.plot_dpi_standard if not dozoom else self.plot_dpi_zoomed + xsize = self.plot_x_size_in * dpi + ysize = xsize // 2 + + self.check_result() event_metadata = self.get_event_metadata() unique_id = f'{str(event_metadata)}_{self.get_nside_string()}' plot_title = f"Run: {event_metadata.run_id} Event {event_metadata.event_id}: Type: {event_metadata.event_type} MJD: {event_metadata.mjd}" - plot_filename = f"{unique_id}.{'plot_zoomed.' if dozoom else ''}pdf" + plot_filename = f"{unique_id}.{'plot_zoomed_legacy.' if dozoom else ''}pdf" print(f"saving plot to {plot_filename}") nsides = self.nsides print(f"available nsides: {nsides}") - maps = [] min_value = np.nan max_value = np.nan minRA=0. @@ -699,12 +702,12 @@ def upload_func(file_buffer, name, title): print(f"preparing plot: {plot_filename}...") # the color map to use - cmap = matplotlib.cm.plasma_r + cmap = self.plot_colormap cmap.set_under(alpha=0.) # make underflows transparent cmap.set_bad(alpha=1., color=(1.,0.,0.)) # make NaNs bright red # prepare the figure canvas - fig = matplotlib.pyplot.figure(figsize=[x_inches,y_inches]) + fig = matplotlib.pyplot.figure(figsize=(self.plot_x_size_in,self.plot_y_size_in)) if dozoom: ax = fig.add_subplot(111) #,projection='cartesian') else: @@ -716,18 +719,7 @@ def upload_func(file_buffer, name, title): image = ax.pcolormesh(ra, dec, grid_map, vmin=min_value, vmax=max_value, rasterized=True, cmap=cmap) # ax.set_xlim(np.pi, -np.pi) - # Use Green's theorem to compute the area - # enclosed by the given contour. - def area(vs): - a = 0 - x0,y0 = vs[0] - for [x1,y1] in vs[1:]: - dx = x1-x0 - dy = y1-y0 - a += 0.5*(y0*dx - x0*dy) - x0 = x1 - y0 = y1 - return a + contour_levels = (np.array([1.39, 4.61, 11.83, 28.74])+min_value)[:2] contour_labels = [r'50%', r'90%', r'3$\sigma$', r'5$\sigma$'][:2] @@ -742,7 +734,9 @@ def area(vs): if not dozoom: # graticule + # mypy error: "Axes" has no attribute "set_longitude_grid" [attr-defined] ax.set_longitude_grid(30) + # mypy error: "Axes" has no attribute "set_latitude_grid" [attr-defined] ax.set_latitude_grid(30) cb = fig.colorbar(image, orientation='horizontal', shrink=.6, pad=0.05, ticks=[min_value, max_value]) cb.ax.xaxis.set_label_text(r"$-2 \ln(L)$") @@ -756,7 +750,8 @@ def area(vs): for i in range(len(contour_labels)): vs = cs_collections[i].get_paths()[0].vertices # Compute area enclosed by vertices. - a = area(vs) # will be in square-radians + # Take absolute values to be independent of orientation of the boundary integral. + a = abs(self.calculate_area(vs)) # will be in square-radians a = a*(180.*180.)/(np.pi*np.pi) # convert to square-degrees leg_labels.append(f'{contour_labels[i]} - area: {a:.2f}sqdeg') @@ -768,6 +763,7 @@ def area(vs): x_width = 1.6 * np.sqrt(a) if np.isnan(x_width): + # error: "QuadContourSet" has no attribute "allsegs" [attr-defined] x_width = 1.6*(max(CS.allsegs[i][0][:,0]) - min(CS.allsegs[i][0][:,0])) print(x_width) y_width = 0.5 * x_width @@ -777,8 +773,8 @@ def area(vs): lower_y = max(minDec -y_width*np.pi/180., -np.pi/2.) upper_y = min(minDec + y_width*np.pi/180., np.pi/2.) - ax.set_xlim( [lower_x, upper_x][::-1]) - ax.set_ylim( [lower_y, upper_y]) + ax.set_xlim(upper_x, lower_x) + ax.set_ylim(lower_y, upper_y) ax.xaxis.set_major_formatter(DecFormatter()) ax.yaxis.set_major_formatter(DecFormatter()) @@ -793,7 +789,10 @@ def area(vs): # cb.ax.xaxis.labelpad = -8 # workaround for issue with viewers, see colorbar docstring - cb.solids.set_edgecolor("face") + # mypy compliance: since cb.solids could be None, we check that it is actually + # a valid object before accessing it + if isinstance(cb.solids, matplotlib.collections.QuadMesh): + cb.solids.set_edgecolor("face") if dozoom: ax.set_aspect('equal') @@ -805,7 +804,9 @@ def area(vs): # Otherwise, add the path effects. effects = [patheffects.withStroke(linewidth=1.1, foreground='w')] + # mypy warnings for artist in ax.findobj(text.Text): + # mypy error: Argument 1 to "set_path_effects" of "Artist" has incompatible type "list[withStroke]"; expected "list[AbstractPathEffect]" [arg-type] artist.set_path_effects(effects) # remove white space around figure @@ -818,44 +819,21 @@ def area(vs): # set the title fig.suptitle(plot_title) - if dosave: - print(f"saving: {plot_filename}...") - - fig.savefig(plot_filename, dpi=dpi, transparent=True) + print(f"saving: {plot_filename}...") - # use io.BytesIO to save this into a memory buffer - imgdata = io.BytesIO() - fig.savefig(imgdata, format='png', dpi=dpi, transparent=True) - imgdata.seek(0) + fig.savefig(plot_filename, dpi=dpi, transparent=True) print("done.") - return imgdata - def create_plot_zoomed(self, - dosave=False, - log_func=None, - upload_func=None, extra_ra=np.nan, extra_dec=np.nan, extra_radius=np.nan, systematics=False, plot_bounding_box=False, - plot_4fgl=False, - final_channels=None): + plot_4fgl=False): """Uses healpy to plot a map.""" - if log_func is None: - def log_func(x): - print(x) - - if upload_func is None: - def upload_func(file_buffer, name, title): - pass - - if final_channels is None: - final_channels=["#test_messaging"] - def bounding_box(ra, dec, theta, phi): shift = ra-180 @@ -865,18 +843,12 @@ def bounding_box(ra, dec, theta, phi): dec_minus = (np.pi/2-np.max(theta))*180./np.pi - dec return ra_plus, ra_minus, dec_plus, dec_minus - y_inches = 3.85 - x_inches = 6. - dpi = 1200. - xsize = x_inches*dpi - ysize = xsize/2. + dpi = self.plot_dpi_zoomed lonra=[-10.,10.] latra=[-10.,10.] - for k in self.result: - if "nside-" not in k: - raise RuntimeError("\"nside\" not in result file..") + self.check_result() event_metadata = self.get_event_metadata() unique_id = f'{str(event_metadata)}_{self.get_nside_string()}' @@ -922,14 +894,14 @@ def bounding_box(ra, dec, theta, phi): grid_map[(tmp_dec, tmp_ra)] = value print("done with map for nside {0}...".format(nside)) - grid_dec = []; grid_ra = []; grid_value = [] + grid_dec_list, grid_ra_list, grid_value_list = [], [], [] for (dec, ra), value in grid_map.items(): - grid_dec.append(dec); grid_ra.append(ra) - grid_value.append(value) - grid_dec = np.asarray(grid_dec) - grid_ra = np.asarray(grid_ra) - grid_value = np.asarray(grid_value) + grid_dec_list.append(dec); grid_ra_list.append(ra) + grid_value_list.append(value) + grid_dec: np.ndarray = np.asarray(grid_dec_list) + grid_ra: np.ndarray = np.asarray(grid_ra_list) + grid_value: np.ndarray = np.asarray(grid_value_list) sorting_indices = np.argsort(grid_value) grid_value = grid_value[sorting_indices] @@ -957,7 +929,7 @@ def bounding_box(ra, dec, theta, phi): print("preparing plot: {0}...".format(plot_filename)) - cmap = matplotlib.cm.plasma_r + cmap = self.plot_colormap cmap.set_under('w') cmap.set_bad(alpha=1., color=(1.,0.,0.)) # make NaNs bright red @@ -1021,34 +993,20 @@ def bounding_box(ra, dec, theta, phi): healpy.projplot(np.pi/2 - minDec, minRA, '*', ms=5, label=r'scan best fit', color='black', zorder=2) - # Use Green's theorem to compute the area - # enclosed by the given contour. - def area(vs): - a = 0 - x0,y0 = vs[0] - for [x1,y1] in vs[1:]: - dx = x1-x0 - dy = y1-y0 - a += 0.5*(y0*dx - x0*dy) - x0 = x1 - y0 = y1 - return a - # Plot the contours contour_areas=[] for contour_level, contour_label, contour_color, contours in zip(contour_levels, contour_labels, contour_colors, contours_by_level): - contour_area = 0 + contour_area = 0. for contour in contours: _ = contour.copy() _[:,1] += np.pi-np.radians(ra) _[:,1] %= 2*np.pi - contour_area += area(_) - contour_area = abs(contour_area) - contour_area *= (180.*180.)/(np.pi*np.pi) # convert to square-degrees - contour_areas.append(contour_area) + contour_area += self.calculate_area(_) + contour_area_sqdeg = abs(contour_area) * (180.*180.)/(np.pi*np.pi) # convert to square-degrees + contour_areas.append(contour_area_sqdeg) contour_label = contour_label + ' - area: {0:.2f} sqdeg'.format( - contour_area) + contour_area_sqdeg) first = True for contour in contours: theta, phi = contour.T @@ -1096,35 +1054,36 @@ def area(vs): ra, ra_plus, np.abs(ra_minus)) + " \n" + \ "\t Dec = {0:.2f} + {1:.2f} - {2:.2f}".format( dec, dec_plus, np.abs(dec_minus)) - log_func(contain_txt) + print(contain_txt) if plot_bounding_box: - bounding_ras = []; bounding_decs = [] + bounding_ras_list, bounding_decs_list = [], [] # lower bound - bounding_ras.extend(list(np.linspace(ra+ra_minus, + bounding_ras_list.extend(list(np.linspace(ra+ra_minus, ra+ra_plus, 10))) - bounding_decs.extend([dec+dec_minus]*10) + bounding_decs_list.extend([dec+dec_minus]*10) # right bound - bounding_ras.extend([ra+ra_plus]*10) - bounding_decs.extend(list(np.linspace(dec+dec_minus, + bounding_ras_list.extend([ra+ra_plus]*10) + bounding_decs_list.extend(list(np.linspace(dec+dec_minus, dec+dec_plus, 10))) # upper bound - bounding_ras.extend(list(np.linspace(ra+ra_plus, + bounding_ras_list.extend(list(np.linspace(ra+ra_plus, ra+ra_minus, 10))) - bounding_decs.extend([dec+dec_plus]*10) + bounding_decs_list.extend([dec+dec_plus]*10) # left bound - bounding_ras.extend([ra+ra_minus]*10) - bounding_decs.extend(list(np.linspace(dec+dec_plus, + bounding_ras_list.extend([ra+ra_minus]*10) + bounding_decs_list.extend(list(np.linspace(dec+dec_plus, dec+dec_minus,10))) # join end to beginning - bounding_ras.append(bounding_ras[0]) - bounding_decs.append(bounding_decs[0]) - bounding_ras = np.asarray(bounding_ras) - bounding_decs = np.asarray(bounding_decs) + bounding_ras_list.append(bounding_ras_list[0]) + bounding_decs_list.append(bounding_decs_list[0]) + + bounding_ras: np.ndarray = np.asarray(bounding_ras_list) + bounding_decs: np.ndarray = np.asarray(bounding_decs_list) bounding_phi = np.radians(bounding_ras) bounding_theta = np.pi/2 - np.radians(bounding_decs) bounding_contour = np.array([bounding_theta, bounding_phi]) bounding_contour_area = 0. - bounding_contour_area = area(bounding_contour.T) + bounding_contour_area = abs(self.calculate_area(bounding_contour.T)) bounding_contour_area *= (180.*180.)/(np.pi*np.pi) # convert to square-degrees contour_label = r'90% Bounding rectangle' + ' - area: {0:.2f} sqdeg'.format( bounding_contour_area) @@ -1132,7 +1091,7 @@ def area(vs): c='r', linestyle='dashed', label=contour_label) # Output contours in RA, dec instead of theta, phi - saving_contours = [] + saving_contours: list = [] for contours in contours_by_level: saving_contours.append([]) for contour in contours: @@ -1150,19 +1109,11 @@ def area(vs): tab = {"ra (rad)": ras, "dec (rad)": decs} savename = unique_id + ".contour_" + val + ".txt" try: - ascii.write(tab, savename, overwrite=True) print("Dumping to", savename) - for i, ch in enumerate(final_channels): - output = io.StringIO() - #output = str.encode(savename) - if dosave: - ascii.write(tab, output, overwrite=True) - output.seek(0) - print(upload_func(output, savename, savename)) - output.truncate(0) - del output - except OSError: - log_func("Memory Error prevented contours from being written") + ascii.write(tab, savename, overwrite=True) + except OSError as err: + print("OS Error prevented contours from being written, maybe a memory issue.") + print(err) uncertainty = [(ra_minus, ra_plus), (dec_minus, dec_plus)] fits_header = format_fits_header( @@ -1235,38 +1186,25 @@ def circular_contour(ra, dec, sigma, nside): print("Contour Area (90%):", contour_areas[1], "degrees (cartesian)", healpy_area, "degrees (scaled)") - if dosave: - # Dump the whole contour - path = unique_id + ".contour.pkl" - print("Saving contour to", path) - with open(path, "wb") as f: - pickle.dump(saving_contours, f) - healpy.write_map(f"{unique_id}.skymap_nside_{mmap_nside}.fits.gz", - equatorial_map, coord = 'C', column_names = ['2LLH'], - extra_header = fits_header, overwrite=True) + # Dump the whole contour + path = unique_id + ".contour.pkl" + print("Saving contour to", path) + with open(path, "wb") as f: + pickle.dump(saving_contours, f) - # Save the figure - print("saving: {0}...".format(plot_filename)) - #ax.invert_xaxis() - fig.savefig(plot_filename, dpi=dpi, transparent=True) + healpy.write_map(f"{unique_id}.skymap_nside_{mmap_nside}.fits.gz", + equatorial_map, coord = 'C', column_names = ['2LLH'], + extra_header = fits_header, overwrite=True) - print("done.") + # Save the figure + print("saving: {0}...".format(plot_filename)) + #ax.invert_xaxis() + fig.savefig(plot_filename, dpi=dpi, transparent=True) - if systematics is True: - title = "Millipede contour, assuming IC160427A systematics:" - else: - title = "Millipede contour, assuming Wilks' Theorem:" - - for i, ch in enumerate(final_channels): - imgdata = io.BytesIO() - fig.savefig(imgdata, format='png', dpi=600, transparent=True) - imgdata.seek(0) + print("done.") - savename = plot_filename[:-4] + ".png" - print(savename) - # config.slack_channel=ch - upload_func(imgdata, savename, title) + savename = plot_filename[:-4] + ".png" + print(savename) plt.close() - return imgdata