From 061f148860a54790372be84c9041d55987d103e1 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Sun, 19 Sep 2021 11:22:48 +0100 Subject: [PATCH 1/8] Add ILOFAR Fido client --- radiospectra/net/__init__.py | 4 +- radiospectra/net/sources/ilofar.py | 85 +++++ .../net/sources/tests/data/ilofar_resp1.html | 305 ++++++++++++++++ .../net/sources/tests/data/ilofar_resp2.html | 337 ++++++++++++++++++ radiospectra/net/sources/tests/test_ilofar.py | 57 +++ 5 files changed, 787 insertions(+), 1 deletion(-) create mode 100644 radiospectra/net/sources/ilofar.py create mode 100644 radiospectra/net/sources/tests/data/ilofar_resp1.html create mode 100644 radiospectra/net/sources/tests/data/ilofar_resp2.html create mode 100644 radiospectra/net/sources/tests/test_ilofar.py diff --git a/radiospectra/net/__init__.py b/radiospectra/net/__init__.py index b868a50..ac86816 100644 --- a/radiospectra/net/__init__.py +++ b/radiospectra/net/__init__.py @@ -1,9 +1,11 @@ from radiospectra.net.attrs import * from radiospectra.net.sources.callisto import CALLISTOClient from radiospectra.net.sources.eovsa import EOVSAClient +from radiospectra.net.sources.ilofar import ILOFARClient from radiospectra.net.sources.psp import RFSClient from radiospectra.net.sources.rstn import RSTNClient from radiospectra.net.sources.stereo import SWAVESClient from radiospectra.net.sources.wind import WAVESClient -__all__ = ["CALLISTOClient", "EOVSAClient", "RFSClient", "SWAVESClient", "RSTNClient", "WAVESClient"] +__all__ = ['CALLISTOClient', 'EOVSAClient', 'RFSClient', 'SWAVESClient', 'RSTNClient', + 'WAVESClient', 'ILOFARClient'] diff --git a/radiospectra/net/sources/ilofar.py b/radiospectra/net/sources/ilofar.py new file mode 100644 index 0000000..3e42f7a --- /dev/null +++ b/radiospectra/net/sources/ilofar.py @@ -0,0 +1,85 @@ +import astropy.units as u +from sunpy.net import attrs as a +from sunpy.net.dataretriever.client import GenericClient, QueryResponse +from sunpy.time import TimeRange +from sunpy.util.scraper import Scraper + +__all__ = ['ILOFARClient'] + + +RECEIVER_FREQUENCIES = { + 'rad1': a.Wavelength(20*u.kHz, 1040*u.kHz), + 'rad2': a.Wavelength(1.075*u.MHz, 13.825*u.MHz), +} + +RECEIVER_EXT = { + 'rad1': 'R1', + 'rad2': 'R2', +} + + +class ILOFARClient(GenericClient): + """ + Provides access to I-LOFAR mode 357 observations from the + data `archive `__ + + Examples + -------- + >>> import radiospectra.net + >>> from sunpy.net import Fido, attrs as a + >>> results = Fido.search(a.Time("2010/10/01", "2010/10/02"), + ... a.Instrument('WAVES')) # doctest: +REMOTE_DATA + >>> results #doctest: +REMOTE_DATA + + Results from 1 Provider: + + 4 Results from the WAVESClient: + Start Time End Time ... Provider Wavelength [2] + ... kHz + ----------------------- ----------------------- ... -------- ----------------- + 2010-10-01 00:00:00.000 2010-10-01 23:59:59.999 ... NASA 20.0 .. 1040.0 + 2010-10-02 00:00:00.000 2010-10-02 23:59:59.999 ... NASA 20.0 .. 1040.0 + 2010-10-01 00:00:00.000 2010-10-01 23:59:59.999 ... NASA 1075.0 .. 13825.0 + 2010-10-02 00:00:00.000 2010-10-02 23:59:59.999 ... NASA 1075.0 .. 13825.0 + + """ + + baseurl = (r'https://data.lofar.ie/%Y/%m/%d/bst/kbt/rcu357_1beam/' + r'%Y%m%d_\d{6}_bst_\d{2}\S{1}.dat') + + pattern = r'{}/{year:4d}{month:2d}{day:2d}_{hour:2d}{minute:2d}{second:2d}' \ + r'_bst_{num:2d}{Polarisation}.dat' + + def search(self, *args, **kwargs): + """ + Query this client for a list of results. + + Parameters + ---------- + *args: `tuple` + `sunpy.net.attrs` objects representing the query. + **kwargs: `dict` + Any extra keywords to refine the search. + + Returns + ------- + A `QueryResponse` instance containing the query result. + """ + matchdict = self._get_match_dict(*args, **kwargs) + metalist = [] + tr = TimeRange(matchdict['Start Time'], matchdict['End Time']) + scraper = Scraper(self.baseurl, regex=True) + filesmeta = scraper._extract_files_meta(tr, extractor=self.pattern) + for i in filesmeta: + rowdict = self.post_search_hook(i, matchdict) + metalist.append(rowdict) + + return QueryResponse(metalist, client=self) + + @classmethod + def register_values(cls): + adict = {a.Instrument: [('ILOFAR', 'Irish LOFAR STATION (IE63)')], + a.Source: [('ILOFAR', 'Irish LOFAR Data Archive')], + a.Provider: [('ILOFAR', 'Irish LOFAR Data Archive')], + a.Wavelength: [('*')]} + return adict diff --git a/radiospectra/net/sources/tests/data/ilofar_resp1.html b/radiospectra/net/sources/tests/data/ilofar_resp1.html new file mode 100644 index 0000000..c1d52ed --- /dev/null +++ b/radiospectra/net/sources/tests/data/ilofar_resp1.html @@ -0,0 +1,305 @@ + diff --git a/radiospectra/net/sources/tests/data/ilofar_resp2.html b/radiospectra/net/sources/tests/data/ilofar_resp2.html new file mode 100644 index 0000000..b86b1dc --- /dev/null +++ b/radiospectra/net/sources/tests/data/ilofar_resp2.html @@ -0,0 +1,337 @@ + diff --git a/radiospectra/net/sources/tests/test_ilofar.py b/radiospectra/net/sources/tests/test_ilofar.py new file mode 100644 index 0000000..cc764f9 --- /dev/null +++ b/radiospectra/net/sources/tests/test_ilofar.py @@ -0,0 +1,57 @@ +from pathlib import Path +from unittest import mock + +import pytest + +from sunpy.net import Fido +from sunpy.net import attrs as a + +from radiospectra.net.sources.ilofar import ILOFARClient + + +@pytest.fixture +def client(): + return ILOFARClient() + + +@pytest.fixture +def html_responses(): + paths = [Path(__file__).parent / 'data' / n for n in ['ilofar_resp1.html', 'ilofar_resp2.html']] + response_htmls = [] + for p in paths: + with p.open('r') as f: + response_htmls.append(f.read()) + + return response_htmls + + +@mock.patch('sunpy.util.scraper.urlopen') +def test_ilofar_client(mock_urlopen, client, html_responses): + mock_urlopen.return_value.read = mock.MagicMock() + mock_urlopen.return_value.read.side_effect = html_responses + mock_urlopen.close = mock.MagicMock(return_value=None) + atr = a.Time('2018/06/01', '2018/06/02') + query = client.search(atr) + + called_urls = ['https://data.lofar.ie/2018/06/01/bst/kbt/rcu357_1beam/', + 'https://data.lofar.ie/2018/06/02/bst/kbt/rcu357_1beam/'] + assert called_urls == [call[0][0] for call in mock_urlopen.call_args_list] + assert len(query) == 8 + assert query[0]['Source'] == 'ILOFAR' + assert query[0]['Provider'] == 'ILOFAR' + assert query[0]['Start Time'].iso == '2018-06-01 10:00:41.000' + assert query[0]['Polarisation'] == 'X' + + +@pytest.mark.remote_data +def test_fido(): + atr = a.Time('2018/06/01', '2018/06/02') + query = Fido.search(atr, a.Instrument('ILOFAR')) + + assert isinstance(query[0].client, ILOFARClient) + query = query[0] + assert len(query) == 8 + assert query[0]['Source'] == 'ILOFAR' + assert query[0]['Provider'] == 'ILOFAR' + assert query[0]['Start Time'].iso == '2018-06-01 10:00:41.000' + assert query[0]['Polarisation'] == 'X' From 48d6292c8e3d1abb2c897104d58feee50a94ec4d Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Mon, 20 Sep 2021 11:40:39 +0100 Subject: [PATCH 2/8] Add ILOFAR spectrogram. --- radiospectra/spectrogram2/__init__.py | 13 + radiospectra/spectrogram2/sources.py | 219 +++++ radiospectra/spectrogram2/spectrogram.py | 764 ++++++++++++++++++ .../spectrogram2/tests/test_ilofar.py | 28 + 4 files changed, 1024 insertions(+) create mode 100644 radiospectra/spectrogram2/__init__.py create mode 100644 radiospectra/spectrogram2/sources.py create mode 100644 radiospectra/spectrogram2/spectrogram.py create mode 100644 radiospectra/spectrogram2/tests/test_ilofar.py diff --git a/radiospectra/spectrogram2/__init__.py b/radiospectra/spectrogram2/__init__.py new file mode 100644 index 0000000..f9bda63 --- /dev/null +++ b/radiospectra/spectrogram2/__init__.py @@ -0,0 +1,13 @@ +""" +The spectrogram2 module aims to provide a more sunpy-like interface to radiospectra data similar to +that of `~sunpy.map.Map` and `~sunpy.timeseries.TimeSeries`. +""" + +from radiospectra.spectrogram2.sources import ( + CALISTOSpectrogram, + EOVSASpectrogram, + ILOFARSpectrogram, + RFSSpectrogram, + SWAVESSpectrogram, +) +from radiospectra.spectrogram2.spectrogram import Spectrogram diff --git a/radiospectra/spectrogram2/sources.py b/radiospectra/spectrogram2/sources.py new file mode 100644 index 0000000..2e3ccea --- /dev/null +++ b/radiospectra/spectrogram2/sources.py @@ -0,0 +1,219 @@ + +import astropy.units as u +from astropy.coordinates.earth import EarthLocation + +from radiospectra.spectrogram2.spectrogram import GenericSpectrogram + +__all__ = ['SWAVESSpectrogram', 'RFSSpectrogram', 'CALISTOSpectrogram', 'EOVSASpectrogram', + 'RSTNSpectrogram', 'ILOFARSpectrogram'] + + +class SWAVESSpectrogram(GenericSpectrogram): + """ + STEREO Waves or S/WAVES, SWAVES Spectrogram + + Examples + -------- + >>> import radiospectra.net + >>> from sunpy.net import Fido, attrs as a + >>> from radiospectra.spectrogram2 import Spectrogram + >>> from radiospectra.net import attrs as ra + >>> query = Fido.search(a.Time('2019/10/05 23:00', '2019/10/06 00:59'), #doctest: +REMOTE_DATA + ... a.Instrument('SWAVES')) #doctest: +REMOTE_DATA + >>> downloaded = Fido.fetch(query[1][0]) #doctest: +REMOTE_DATA + >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA + >>> spec #doctest: +REMOTE_DATA + + >>> spec.plot() #doctest: +REMOTE_DATA + """ + def __init__(self, data, meta, **kwargs): + super().__init__(meta=meta, data=data, **kwargs) + + @property + def receiver(self): + """ + The name of the receiver + """ + return self.meta['receiver'] + + @classmethod + def is_datasource_for(cls, data, meta, **kwargs): + return meta['instrument'] == 'swaves' + + +class WAVESSpectrogram(GenericSpectrogram): + """ + Wind Waves Spectrogram + + Examples + -------- + >>> import radiospectra.net + >>> from sunpy.net import Fido, attrs as a + >>> from radiospectra.spectrogram2 import Spectrogram + >>> from radiospectra.net import attrs as ra + >>> query = Fido.search(a.Time('2019/10/05 23:00', '2019/10/06 00:59'), #doctest: +REMOTE_DATA + ... a.Instrument('WAVES')) #doctest: +REMOTE_DATA + >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA + >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA + >>> spec #doctest: +REMOTE_DATA + + >>> spec.plot() #doctest: +REMOTE_DATA + """ + def __init__(self, data, meta, **kwargs): + super().__init__(meta=meta, data=data, **kwargs) + + @property + def receiver(self): + """ + The name of the receiver + """ + return self.meta['receiver'] + + @property + def background(self): + """ + The background subtracted from the data + """ + return self.meta.bg + + @classmethod + def is_datasource_for(cls, data, meta, **kwargs): + return meta.get('instrument', None) == 'WAVES' + + +class RFSSpectrogram(GenericSpectrogram): + """ + Parker Solar Probe FIELDS/Radio Frequency Spectrometer (RFS) Spectrogram + + >>> import radiospectra.net + >>> from sunpy.net import Fido, attrs as a + >>> from radiospectra.spectrogram2 import Spectrogram + >>> from radiospectra.net import attrs as ra + >>> query = Fido.search(a.Time('2019/10/05 23:00', '2019/10/06 00:59'), #doctest: +REMOTE_DATA + ... a.Instrument.rfs) #doctest: +REMOTE_DATA + >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA + >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA + >>> spec #doctest: +REMOTE_DATA + + >>> spec.plot() #doctest: +REMOTE_DATA + """ + def __init__(self, data, meta, **kwargs): + super().__init__(meta=meta, data=data, **kwargs) + + @property + def level(self): + return self.meta['cdf_meta']['Data_type'].split('>')[0] + + @property + def version(self): + return int(self.meta['cdf_meta']['Data_version']) + + @classmethod + def is_datasource_for(cls, data, meta, **kwargs): + return (meta['observatory'] == 'PSP' and meta['instrument'] == 'FIELDS/RFS' + and meta['detector'] in ('lfr', 'hfr')) + + +class CALISTOSpectrogram(GenericSpectrogram): + """ + CALISTO Spectrogram from the e-CALISTO network + + Examples + -------- + >>> import radiospectra.net + >>> from sunpy.net import Fido, attrs as a + >>> from radiospectra.spectrogram2 import Spectrogram + >>> from radiospectra.net import attrs as ra + >>> query = Fido.search(a.Time('2019/10/05 23:00', '2019/10/06 00:59'), #doctest: +REMOTE_DATA + ... a.Instrument('eCALLISTO'), ra.Observatory('ALASKA')) #doctest: +REMOTE_DATA + >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA + >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA + >>> spec #doctest: +REMOTE_DATA + + >>> spec.plot() #doctest: +REMOTE_DATA + """ + def __init__(self, data, meta, **kwargs): + super().__init__(meta=meta, data=data, **kwargs) + + @property + def observatory_location(self): + lat = self.meta['fits_meta']['OBS_LAT'] * u.deg * 1.0 if \ + self.meta['fits_meta']['OBS_LAC'] == 'N' else -1.0 + lon = self.meta['fits_meta']['OBS_LON'] * u.deg * 1.0 if \ + self.meta['fits_meta']['OBS_LOC'] == 'E' else -1.0 + height = self.meta['fits_meta']['OBS_ALT'] * u.m + return EarthLocation(lat=lat, lon=lon, height=height) + + @classmethod + def is_datasource_for(cls, data, meta, **kwargs): + return meta['instrument'] == 'e-CALLISTO' or meta['detector'] == 'e-CALLISTO' + + +class EOVSASpectrogram(GenericSpectrogram): + """ + Extend Owen Valley Array (EOVSA) Spectrogram + + Examples + -------- + >>> import radiospectra.net + >>> from sunpy.net import Fido, attrs as a + >>> from radiospectra.spectrogram2 import Spectrogram + >>> query = Fido.search(a.Time('2021/05/07 00:00', '2021/05/07 23:00'), a.Instrument.eovsa) #doctest: +REMOTE_DATA + >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA + >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA + >>> spec #doctest: +REMOTE_DATA + + >>> spec.plot() #doctest: +REMOTE_DATA + """ + def __init__(self, data, meta, **kwargs): + super().__init__(meta=meta, data=data, **kwargs) + + @property + def polarisation(self): + return self.meta['fits_meta']['POLARIZA'] + + @classmethod + def is_datasource_for(cls, data, meta, **kwargs): + return meta['instrument'] == 'EOVSA' or meta['detector'] == 'EOVSA' + + # TODO fix time gaps for plots need to render them as gaps + # can prob do when generateing proper pcolormesh gird but then prob doesn't belong here + + +class RSTNSpectrogram(GenericSpectrogram): + """ + Radio Solar Telescope Network + + Examples + -------- + >>> import radiospectra.net + >>> from sunpy.net import Fido, attrs as a + >>> from radiospectra.spectrogram2 import Spectrogram + >>> query = Fido.search(a.Time('2017/09/07 00:00', '2017/09/07 23:00'), a.Instrument.rstn) #doctest: +REMOTE_DATA + >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA + >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA + >>> spec #doctest: +REMOTE_DATA + + >>> spec.plot() #doctest: +REMOTE_DATA + """ + + @classmethod + def is_datasource_for(cls, data, meta, **kwargs): + return meta['instrument'] == 'RSTN' + + +class ILOFARSpectrogram(GenericSpectrogram): + """ + Irish LOFAR Station + """ + @property + def mode(self): + return self.meta.get('mode') + + @property + def polarisation(self): + return self.meta.get('polarisation') + + @classmethod + def is_datasource_for(cls, data, meta, **kwargs): + return meta['instrument'] == 'ILOFAR' diff --git a/radiospectra/spectrogram2/spectrogram.py b/radiospectra/spectrogram2/spectrogram.py new file mode 100644 index 0000000..1ead98a --- /dev/null +++ b/radiospectra/spectrogram2/spectrogram.py @@ -0,0 +1,764 @@ +import os +import glob +import gzip +import struct +import pathlib +import warnings +from pathlib import Path +from collections import OrderedDict +from urllib.request import Request, urlopen + +import cdflib +import matplotlib.dates as mdates +import numpy as np +import pandas as pd +from matplotlib import pyplot as plt +from matplotlib.image import NonUniformImage +from scipy.io.idl import readsav + +import astropy.units as u +from astropy.io.fits import Header +from astropy.time import Time +from astropy.visualization import quantity_support +from sunpy.data import cache +from sunpy.io import fits +from sunpy.net import attrs as a +from sunpy.time import parse_time +from sunpy.util.datatype_factory_base import ( + BasicRegistrationFactory, + MultipleMatchError, + NoMatchError, + ValidationFunctionError, +) +from sunpy.util.exceptions import SunpyUserWarning +from sunpy.util.functools import seconddispatch +from sunpy.util.metadata import MetaDict +from sunpy.util.util import expand_list + +SUPPORTED_ARRAY_TYPES = (np.ndarray,) +try: + import dask.array + SUPPORTED_ARRAY_TYPES += (dask.array.Array,) +except ImportError: + pass + +quantity_support() + +__all__ = ['SpectrogramFactory', 'Spectrogram', 'GenericSpectrogram'] + + +class NoSpectrogramInFileError(Exception): + pass + + +class SpectraMetaValidationError(AttributeError): + pass + + +# Can use sunpy.until.io once min version has these +def is_url(arg): + try: + urlopen(arg) + except Exception: + return False + return True + + +# In python<3.8 paths with un-representable chars (ie. '*' on windows) +# raise an error, so make our own version that returns False instead of +# erroring. These can be removed when we support python >= 3.8 +# https://docs.python.org/3/library/pathlib.html#methods +def is_file(path): + try: + return path.is_file() + except Exception: + return False + + +def is_dir(path): + try: + return path.is_dir() + except Exception: + return False + + +def possibly_a_path(obj): + """ + Check if ``obj`` can be coerced into a pathlib.Path object. + Does *not* check if the path exists. + """ + try: + pathlib.Path(obj) + return True + except Exception: + return False + + +def parse_path(path, f, **kwargs): + """ + Read in a series of files at *path* using the function *f*. + + Parameters + ---------- + path : pathlib.Path + f : callable + Must return a list of read-in data. + kwargs : + Additional keyword arguments are handed to ``f``. + + Returns + ------- + list + List of files read in by ``f``. + """ + if not isinstance(path, os.PathLike): + raise ValueError('path must be a pathlib.Path object') + path = path.expanduser() + if is_file(path): + return list(f(path, **kwargs)) + elif is_dir(path): + read_files = [] + for afile in sorted(path.glob('*')): + read_files += f(afile, **kwargs) + return read_files + elif glob.glob(str(path)): + read_files = [] + for afile in sorted(glob.glob(str(path))): + read_files += f(afile, **kwargs) + return read_files + else: + raise ValueError(f'Did not find any files at {path}') + + +class PcolormeshPlotMixin: + """ + Class provides plotting functions using `~pcolormesh`. + """ + def plot(self, axes=None, **kwargs): + """ + Plot the spectrogram. + + Parameters + ---------- + axes : `matplotlib.axis.Axes`, optional + The axes where the plot will be added. + kwargs : + Arguments pass to the plot call `pcolormesh`. + + Returns + ------- + """ + if axes is None: + fig, axes = plt.subplots() + else: + fig = axes.get_figure() + + if hasattr(self.data, 'value'): + data = self.data.value + else: + data = self.data + + title = f'{self.observatory}, {self.instrument}' + if self.instrument != self.detector: + title = f'{title}, {self.detector}' + + axes.set_title(title) + axes.plot(self.times.datetime[[0, -1]], self.frequencies[[0, -1]], + linestyle='None', marker='None') + if (self.times.shape[0] == self.data.shape[0] + and self.frequencies.shape[0] == self.data.shape[1]): + axes.pcolormesh(self.times.datetime, self.frequencies.value, data, shading='auto', + **kwargs) + else: + axes.pcolormesh(self.times.datetime, self.frequencies.value, data[:-1, :-1], + shading='auto', **kwargs) + axes.set_xlim(self.times.datetime[0], self.times.datetime[-1]) + locator = mdates.AutoDateLocator(minticks=4, maxticks=8) + formatter = mdates.ConciseDateFormatter(locator) + axes.xaxis.set_major_locator(locator) + axes.xaxis.set_major_formatter(formatter) + fig.autofmt_xdate() + + +class NonUniformImagePlotMixin: + """ + Class provides plotting functions using `NonUniformImage` + """ + def plotim(self, fig=None, axes=None, **kwargs): + if axes is None: + fig, axes = plt.subplots() + + im = NonUniformImage(axes, interpolation='none', **kwargs) + im.set_data(mdates.date2num(self.times.datetime), self.frequencies.value, self.data) + axes.images.append(im) + + +class GenericSpectrogram(PcolormeshPlotMixin, NonUniformImagePlotMixin): + """ + Base spectrogram class all new spectrograms will inherit. + + Attributes + ---------- + meta : `dict-like` + Meta data for the spectrogram. + data : `numpy.ndarray` + The spectrogram data itself a 2D array. + """ + _registry = {} + + def __init_subclass__(cls, **kwargs): + super().__init_subclass__(**kwargs) + if hasattr(cls, 'is_datasource_for'): + cls._registry[cls] = cls.is_datasource_for + + def __init__(self, data, meta, **kwargs): + self.data = data + self.meta = meta + + self._validate_meta() + + @property + def observatory(self): + """ + The name of the observatory which recorded the spectrogram. + """ + return self.meta['observatory'].upper() + + @property + def instrument(self): + """ + The name of the instrument which recorded the spectrogram. + """ + return self.meta['instrument'].upper() + + @property + def detector(self): + """ + The detector which recorded the spectrogram. + """ + return self.meta['detector'].upper() + + @property + def start_time(self): + """ + The start time of the spectrogram. + """ + return self.meta['start_time'] + + @property + def end_time(self): + """ + The end time of the spectrogram. + """ + return self.meta['end_time'] + + @property + def wavelength(self): + """ + The wavelength range of the spectrogram. + """ + return self.meta['wavelength'] + + @property + def times(self): + """ + The times of the spectrogram. + """ + return self.meta['times'] + + @property + def frequencies(self): + """ + The frequencies of the spectrogram. + """ + return self.meta['freqs'] + + def _validate_meta(self): + """ + Validates the meta-information associated with a Spectrogram. + + This method includes very basic validation checks which apply to + all of the kinds of files that radiospectra can read. Datasource-specific + validation should be handled in the relevant file in the + radiospectra.spectrogram2.sources file. + """ + msg = 'Spectrogram coordinate units for {} axis not present in metadata.' + err_message = [] + for i, ax in enumerate(['times', 'freqs']): + if self.meta.get(ax) is None: + err_message.append(msg.format(ax)) + + if err_message: + raise SpectraMetaValidationError('\n'.join(err_message)) + + def __repr__(self): + return (f'<{self.__class__.__name__} {self.observatory}, {self.instrument}, {self.detector}' + f' {self.wavelength.min} - {self.wavelength.max},' + f' {self.start_time.isot} to {self.end_time.isot}>') + + +class SpectrogramFactory(BasicRegistrationFactory): + """ + A factory for generating spectrograms. + + Parameters + ---------- + \\*inputs + `str` or `pathlib.Path` to the file. + + Returns + ------- + `radiospectra.spectrogram2.Spectrogram` + The spectrogram for the give file + """ + def _validate_meta(self, meta): + """ + Validate a meta argument. + """ + if isinstance(meta, Header): + return True + elif isinstance(meta, dict): + return True + else: + return False + + def _parse_args(self, *args, silence_errors=False, **kwargs): + """ + Parses an args list into data-header pairs. + args can contain any mixture of the following entries: + * tuples of data,header + * data, header not in a tuple + * data, wcs object in a tuple + * data, wcs object not in a tuple + * filename, as a str or pathlib.Path, which will be read + * directory, as a str or pathlib.Path, from which all files will be read + * glob, from which all files will be read + * url, which will be downloaded and read + * lists containing any of the above. + + Example + ------- + self._parse_args(data, header, + (data, header), + ['file1', 'file2', 'file3'], + 'file4', + 'directory1', + '*.fits') + """ + # Account for nested lists of items + args = expand_list(args) + + # Sanitise the input so that each 'type' of input corresponds to a different + # class, so single dispatch can be used later + nargs = len(args) + i = 0 + while i < nargs: + arg = args[i] + if isinstance(arg, SUPPORTED_ARRAY_TYPES): + # The next two items are data and a header + data = args.pop(i) + header = args.pop(i) + args.insert(i, (data, header)) + nargs -= 1 + elif isinstance(arg, str) and is_url(arg): + # Replace URL string with a Request object to dispatch on later + args[i] = Request(arg) + elif possibly_a_path(arg): + # Replace path strings with Path objects + args[i] = pathlib.Path(arg) + i += 1 + + # Parse the arguments + # Note that this list can also contain GenericSpectrogram if they are directly given to + # the factory + data_header_pairs = [] + for arg in args: + try: + data_header_pairs += self._parse_arg(arg, **kwargs) + except NoSpectrogramInFileError as e: + if not silence_errors: + raise + warnings.warn(f"One of the arguments failed to parse with error: {e}", + SunpyUserWarning) + + return data_header_pairs + + # Note that post python 3.8 this can be @functools.singledispatchmethod + @seconddispatch + def _parse_arg(self, arg, **kwargs): + """ + Take a factory input and parse into (data, header) pairs. + Must return a list, even if only one pair is returned. + """ + raise ValueError(f"Invalid input: {arg}") + + @_parse_arg.register(tuple) + def _parse_tuple(self, arg, **kwargs): + # Data-header + data, header = arg + pair = data, header + if self._validate_meta(header): + pair = (data, OrderedDict(header)) + return [pair] + + @_parse_arg.register(GenericSpectrogram) + def _parse_map(self, arg, **kwargs): + return [arg] + + @_parse_arg.register(Request) + def _parse_url(self, arg, **kwargs): + url = arg.full_url + path = str(cache.download(url).absolute()) + pairs = self._read_file(path, **kwargs) + return pairs + + @_parse_arg.register(pathlib.Path) + def _parse_path(self, arg, **kwargs): + return parse_path(arg, self._read_file, **kwargs) + + def __call__(self, *args, silence_errors=False, **kwargs): + """ + Method for running the factory. Takes arbitrary arguments and + keyword arguments and passes them to a sequence of pre-registered types + to determine which is the correct spectrogram-type to build. + Arguments args and kwargs are passed through to the validation + function and to the constructor for the final type. For spectrogram types, + validation function must take a data-header pair as an argument. + + Parameters + ---------- + silence_errors : `bool`, optional + If set, ignore data-header pairs which cause an exception. + Default is ``False``. + + Notes + ----- + Extra keyword arguments are passed through to `sunpy.io.read_file` such + as `memmap` for FITS files. + """ + data_header_pairs = self._parse_args(*args, silence_errors=silence_errors, **kwargs) + new_maps = list() + + # Loop over each registered type and check to see if WidgetType + # matches the arguments. If it does, use that type. + for pair in data_header_pairs: + if isinstance(pair, GenericSpectrogram): + new_maps.append(pair) + continue + data, header = pair + meta = MetaDict(header) + + try: + new_map = self._check_registered_widgets(data, meta, **kwargs) + new_maps.append(new_map) + except (NoMatchError, MultipleMatchError, + ValidationFunctionError, SpectraMetaValidationError) as e: + if not silence_errors: + raise + warnings.warn(f'One of the data, header pairs failed to validate with: {e}', + SunpyUserWarning) + + if not len(new_maps): + raise RuntimeError('No maps loaded') + + if len(new_maps) == 1: + return new_maps[0] + + return new_maps + + def _check_registered_widgets(self, data, meta, **kwargs): + + candidate_widget_types = list() + + for key in self.registry: + + # Call the registered validation function for each registered class + if self.registry[key](data, meta, **kwargs): + candidate_widget_types.append(key) + + n_matches = len(candidate_widget_types) + + if n_matches == 0: + if self.default_widget_type is None: + raise NoMatchError("No types match specified arguments and no default is set.") + else: + candidate_widget_types = [self.default_widget_type] + elif n_matches > 1: + raise MultipleMatchError("Too many candidate types identified " + f"({candidate_widget_types}). " + "Specify enough keywords to guarantee unique type " + "identification.") + + # Only one is found + WidgetType = candidate_widget_types[0] + + return WidgetType(data, meta, **kwargs) + + def _read_file(self, file, **kwargs): + file = Path(file) + extensions = file.suffixes + first_extension = extensions[0].lower() + if first_extension == '.dat': + return self._read_dat(file) + elif first_extension in ('.r1', '.r2'): + return self._read_idl_sav(file, instrument='waves') + elif first_extension == '.cdf': + return self._read_cdf(file) + elif first_extension == '.srs': + return self._read_srs(file) + elif first_extension in ('.fits', '.fit', '.fts', 'fit.gz'): + return self._read_fits(file) + else: + raise ValueError(f'Extension {extensions[0]} not supported.') + + @staticmethod + def _read_dat(file): + if 'swaves' in file.name: + name, prod, date, spacecraft, receiver = file.stem.split('_') + # frequency range + freqs = np.genfromtxt(file, max_rows=1) * u.kHz + # bg which is already subtracted from data + bg = np.genfromtxt(file, skip_header=1, max_rows=1) + # data + data = np.genfromtxt(file, skip_header=2) + times = data[:, 0] * u.min + data = data[:, 1:].T + + meta = {'instrument': name, 'observatory': f'STEREO {spacecraft.upper()}', + 'product': prod, 'start_time': Time.strptime(date, '%Y%m%d'), + 'wavelength': a.Wavelength(freqs[0], freqs[-1]), 'detector': receiver, + 'freqs': freqs, 'background': bg} + + meta['times'] = meta['start_time'] + times + meta['end_time'] = meta['start_time'] + times[-1] + + return data, meta + elif 'bst' in file.name: + def subband_to_freq(subband, obs_mode): + """ + Converts LOFAR single station subbands to frequency + + Parameters + ---------- + subband : `int` + Subband number + obs_mode : `int` + Observation mode 3, 5, 7 + + Return + ------ + freq + frequency as astropy.units.Quantity (MHz) + """ + nyq_zone_dict = {3: 1, 5: 2, 7: 3} + nyq_zone = nyq_zone_dict[obs_mode] + clock_dict = {3: 200, 4: 160, 5: 200, 6: 160, 7: 200} # MHz + clock = clock_dict[obs_mode] + freq = (nyq_zone - 1 + subband / 512) * (clock / 2) + return freq * u.MHz # MHz + + subbands = (np.arange(54, 454, 2), np.arange(54, 454, 2), np.arange(54, 230, 2)) + num_subbands = 488 + + data = np.fromfile(file) + polarisation = file.stem[-1] + + num_times = data.shape[0] / num_subbands + data = data.reshape(-1, num_subbands).T # (Freq x Time).T = (Time x Freq) + dt = np.arange(num_times) * 1 * u.s + start_time = Time.strptime(file.name.split('_bst')[0], "%Y%m%d_%H%M%S") + times = start_time + dt + + obs_mode = (3, 5, 7) + + freqs = [subband_to_freq(sb, mode) for sb, mode in zip(subbands, obs_mode)] + + # 1st 200 sbs mode 3, next 200 sbs mode 5, last 88 sbs mode 7 + spec = {0: data[:200, :], 1: data[200:400, :], 2: data[400:, :]} + data_header_pairs = [] + for i in range(3): + meta = {'instrument': 'ILOFAR', 'observatory': 'Birr (IE613)', + 'start_time': times[0], 'mode': obs_mode[i], + 'wavelength': a.Wavelength(freqs[i][0], freqs[i][-1]), + 'freqs': freqs[i], 'times': times, 'end_time': times[-1], + 'detector': 'ILOFAR', 'polarisation': polarisation} + + data_header_pairs.append((spec[i], meta)) + + return data_header_pairs + + @staticmethod + def _read_srs(file): + + with file.open('rb') as buff: + data = buff.read() + if file.suffixes[-1] == '.gz': + data = gzip.decompress(data) + + # Data is store as a series of records made of different numbers of bytes + # General header information + # 1 Year (last 2 digits) Byte integer (unsigned) + # 2 Month number (1 to 12) " + # 3 Day (1 to 31) " + # 4 Hour (0 to 23 UT) " + # 5 Minute (0 to 59) " + # 6 Second at start of scan (0 to 59) " + # 7 Site Number (0 to 255) " + # 8 Number of bands in the record (2) " + # + # Band 1 (A-band) header information + # 9,10 Start Frequency (MHz) Word integer (16 bits) + # 11,12 End Frequency (MHz) " + # 13,14 Number of bytes in data record (401)" + # 15 Analyser reference level Byte integer + # 16 Analyser attenuation (dB) " + # + # Band 2 (B-band) header information + # 17-24 As for band 1 + # + # Spectrum Analyser data + # 25-425 401 data bytes for band 1 (A-band) + # 426-826 401 data bytes for band 2 (B-band) + record_struc = struct.Struct('B'*8 + 'H'*3 + 'B'*2 + + 'H'*3 + 'B'*2 + 'B' * 401 + 'B' * 401) + + records = record_struc.iter_unpack(data) + + # Map of numeric records to locations + site_map = {1: 'Palehua', 2: 'Holloman', 3: 'Learmonth', 4: 'San Vito'} + + df = pd.DataFrame([(*r[:18], np.array(r[18:419]), np.array(r[419:820])) + for r in records]) + df.columns = ['year', 'month', 'day', 'hour', 'minute', 'second', 'site', ' num_bands', + 'start_freq1', 'end_freq1', 'num_bytes1', 'analyser_ref1', + 'analyser_atten1', 'start_freq2', 'end_freq2', 'num_bytes2', + 'analyser_ref2', 'analyser_atten2', 'spec1', 'spec2'] + + # Hack to make to_datetime work - earliest dates seem to be 2000 and won't be + # around in 3000! + df['year'] = df['year'] + 2000 + df['time'] = pd.to_datetime(df[['year', 'month', 'day', 'hour', 'minute', 'second']]) + + # Equations taken from document + n = np.arange(1, 402) + freq_a = (25 + 50 * (n - 1) / 400) * u.MHz + freq_b = (75 + 105 * (n - 1) / 400) * u.MHz + freqs = np.hstack([freq_a, freq_b]) + + data = np.hstack([np.vstack(df[name].to_numpy()) for name in ['spec1', 'spec2']]).T + times = Time(df['time']) + + meta = {'instrument': 'RSTN', 'observatory': site_map[df['site'][0]], + 'start_time': times[0], 'end_time': times[-1], 'detector': 'RSTN', + 'wavelength': a.Wavelength(freqs[0], freqs[-1]), 'freqs': freqs, 'times': times} + + return data, meta + + @staticmethod + def _read_cdf(file): + cdf = cdflib.CDF(file) + cdf_meta = cdf.globalattsget() + if (cdf_meta.get('Project', '') == 'PSP' + and cdf_meta.get('Source_name') == 'PSP_FLD>Parker Solar Probe FIELDS' + and 'Radio Frequency Spectrometer' in cdf_meta.get('Descriptor')): + short, _long = cdf_meta['Descriptor'].split('>') + detector = short[4:].lower() + + times, data, freqs = [cdf.varget(name) for name in + [f'epoch_{detector}_auto_averages_ch0_V1V2', + f'psp_fld_l2_rfs_{detector}_auto_averages_ch0_V1V2', + f'frequency_{detector}_auto_averages_ch0_V1V2']] + + times = Time(times << u.ns, format='cdf_tt2000') + freqs = freqs[0, :] << u.Hz + data = data.T << u.Unit('Volt**2/Hz') + + meta = { + 'cdf_meta': cdf_meta, + 'detector': detector, + 'instrument': 'FIELDS/RFS', + 'observatory': 'PSP', + 'start_time': times[0], + 'end_time': times[-1], + 'wavelength': a.Wavelength(freqs.min(), freqs.max()), + 'times': times, + 'freqs': freqs + } + return data, meta + + @staticmethod + def _read_fits(file): + hd_pairs = fits.read(file) + + if 'e-CALLISTO' in hd_pairs[0].header.get('CONTENT', ''): + data = hd_pairs[0].data + times = hd_pairs[1].data['TIME'].flatten() * u.s + freqs = hd_pairs[1].data['FREQUENCY'].flatten() * u.MHz + start_time = parse_time(hd_pairs[0].header['DATE-OBS'] + + ' ' + hd_pairs[0].header['TIME-OBS']) + end_time = parse_time(hd_pairs[0].header['DATE-END'] + + ' ' + hd_pairs[0].header['TIME-END']) + times = start_time + times + meta = { + 'fits_meta': hd_pairs[0].header, + 'detector': 'e-CALLISTO', + 'instrument': 'e-CALLISTO', + 'observatory': hd_pairs[0].header['INSTRUME'], + 'start_time': start_time, + 'end_time': end_time, + 'wavelength': a.Wavelength(freqs.min(), freqs.max()), + 'times': times, + 'freqs': freqs + } + return data, meta + elif hd_pairs[0].header.get('TELESCOP', '') == 'EOVSA': + times = Time(hd_pairs[2].data['mjd'] + hd_pairs[2].data['time'] / 1000.0 / 86400., + format='mjd') + freqs = hd_pairs[1].data['sfreq'] * u.GHz + data = hd_pairs[0].data + start_time = parse_time(hd_pairs[0].header['DATE_OBS']) + end_time = parse_time(hd_pairs[0].header['DATE_END']) + + meta = { + 'fits_meta': hd_pairs[0].header, + 'detector': 'EOVSA', + 'instrument': 'EOVSA', + 'observatory': 'Owens Valley', + 'start_time': start_time, + 'end_time': end_time, + 'wavelength': a.Wavelength(freqs.min(), freqs.max()), + 'times': times, + 'freqs': freqs + } + return data, meta + + @staticmethod + def _read_idl_sav(file, instrument=None): + data = readsav(file) + if instrument == 'waves': + # See https://solar-radio.gsfc.nasa.gov/wind/one_minute_doc.html + data_array = data['arrayb'] + # frequency range + if file.suffix == '.R1': + freqs = np.linspace(20, 1040, 256) * u.kHz + receiver = 'RAD1' + elif file.suffix == '.R2': + freqs = np.linspace(1.075, 13.825, 256) * u.MHz + receiver = 'RAD2' + # bg which is already subtracted from data ? + bg = data_array[:, -1] + data = data_array[:, :-1] + + start_time = Time.strptime(file.stem, '%Y%m%d') + end_time = start_time + 86399*u.s + times = start_time + (np.arange(1440) * 60 + 30)*u.s + + meta = {'instrument': 'WAVES', 'observatory': 'WIND', 'start_time': start_time, + 'end_time': end_time, 'wavelength': a.Wavelength(freqs[0], freqs[-1]), + 'detector': receiver, 'freqs': freqs, 'times': times, 'background': bg} + + return data, meta + + +Spectrogram = SpectrogramFactory(registry=GenericSpectrogram._registry, + default_widget_type=GenericSpectrogram) diff --git a/radiospectra/spectrogram2/tests/test_ilofar.py b/radiospectra/spectrogram2/tests/test_ilofar.py new file mode 100644 index 0000000..dd1fa0e --- /dev/null +++ b/radiospectra/spectrogram2/tests/test_ilofar.py @@ -0,0 +1,28 @@ +from pathlib import Path +from unittest import mock + +import numpy as np + +from radiospectra.spectrogram2.spectrogram import Spectrogram + + +@mock.patch('numpy.fromfile') +@mock.patch('radiospectra.spectrogram2.spectrogram.is_file') +def test_ilofar(mock_is_file, mock_fromfile): + mock_fromfile.return_value = np.ones(10117216) + mock_is_file.return_value = True + + spec = Spectrogram(Path('20180602_063247_bst_00X.dat')) + assert len(spec) == 3 + assert spec[0].polarisation == 'X' + assert spec[0].start_time.iso == '2018-06-02 06:32:47.000' + assert spec[0].end_time.iso == '2018-06-02 12:18:18.000' + assert spec[0].mode == 3 + assert spec[0].frequencies[0].to_value('MHz') == 10.546875 + assert spec[0].frequencies[-1].to_value('MHz') == 88.28125 + assert spec[1].mode == 5 + assert spec[1].frequencies[0].to_value('MHz') == 110.546875 + assert spec[1].frequencies[-1].to_value('MHz') == 188.28125 + assert spec[2].mode == 7 + assert spec[2].frequencies[0].to_value('MHz') == 210.546875 + assert spec[2].frequencies[-1].to_value('MHz') == 244.53125 From 565b1f03a3fca41cc8fd23af378bdb1af8fef21e Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Mon, 20 Sep 2021 12:52:00 +0100 Subject: [PATCH 3/8] Rename client and spectrogram add more tests --- changelog/57.feature.rst | 1 + radiospectra/net/__init__.py | 4 +- radiospectra/net/sources/ilofar.py | 53 ++++++++++++++----- radiospectra/net/sources/tests/test_ilofar.py | 45 ++++++++++++++-- radiospectra/spectrogram2/__init__.py | 2 +- radiospectra/spectrogram2/sources.py | 6 +-- 6 files changed, 88 insertions(+), 23 deletions(-) create mode 100644 changelog/57.feature.rst diff --git a/changelog/57.feature.rst b/changelog/57.feature.rst new file mode 100644 index 0000000..58b9927 --- /dev/null +++ b/changelog/57.feature.rst @@ -0,0 +1 @@ +Add `sunpy.net.Fido` client `~radiospectra.net.sources.ilofar.ILOFARMode357` and spectrogram class `~radiospectra.spectrogram2.sources.ILOFARMode357` for ILOFAR mode 357 observations. diff --git a/radiospectra/net/__init__.py b/radiospectra/net/__init__.py index ac86816..1206211 100644 --- a/radiospectra/net/__init__.py +++ b/radiospectra/net/__init__.py @@ -1,11 +1,11 @@ from radiospectra.net.attrs import * from radiospectra.net.sources.callisto import CALLISTOClient from radiospectra.net.sources.eovsa import EOVSAClient -from radiospectra.net.sources.ilofar import ILOFARClient +from radiospectra.net.sources.ilofar import ILOFARMode357Client from radiospectra.net.sources.psp import RFSClient from radiospectra.net.sources.rstn import RSTNClient from radiospectra.net.sources.stereo import SWAVESClient from radiospectra.net.sources.wind import WAVESClient __all__ = ['CALLISTOClient', 'EOVSAClient', 'RFSClient', 'SWAVESClient', 'RSTNClient', - 'WAVESClient', 'ILOFARClient'] + 'WAVESClient', 'ILOFARMode357Client'] diff --git a/radiospectra/net/sources/ilofar.py b/radiospectra/net/sources/ilofar.py index 3e42f7a..772932a 100644 --- a/radiospectra/net/sources/ilofar.py +++ b/radiospectra/net/sources/ilofar.py @@ -1,24 +1,19 @@ +import numpy as np + import astropy.units as u from sunpy.net import attrs as a from sunpy.net.dataretriever.client import GenericClient, QueryResponse from sunpy.time import TimeRange from sunpy.util.scraper import Scraper -__all__ = ['ILOFARClient'] - +__all__ = ['ILOFARMode357Client'] -RECEIVER_FREQUENCIES = { - 'rad1': a.Wavelength(20*u.kHz, 1040*u.kHz), - 'rad2': a.Wavelength(1.075*u.MHz, 13.825*u.MHz), -} +from radiospectra.net.attrs import PolType -RECEIVER_EXT = { - 'rad1': 'R1', - 'rad2': 'R2', -} +RECEIVER_FREQUENCIES = a.Wavelength(10.546875*u.MHz, 244.53125*u.MHz) -class ILOFARClient(GenericClient): +class ILOFARMode357Client(GenericClient): """ Provides access to I-LOFAR mode 357 observations from the data `archive `__ @@ -28,7 +23,7 @@ class ILOFARClient(GenericClient): >>> import radiospectra.net >>> from sunpy.net import Fido, attrs as a >>> results = Fido.search(a.Time("2010/10/01", "2010/10/02"), - ... a.Instrument('WAVES')) # doctest: +REMOTE_DATA + ... a.Instrument('ILOFAR')) # doctest: +REMOTE_DATA >>> results #doctest: +REMOTE_DATA Results from 1 Provider: @@ -50,6 +45,23 @@ class ILOFARClient(GenericClient): pattern = r'{}/{year:4d}{month:2d}{day:2d}_{hour:2d}{minute:2d}{second:2d}' \ r'_bst_{num:2d}{Polarisation}.dat' + @classmethod + def _check_wavelengths(cls, wavelength): + """ + Check for overlap between given wavelength and receiver frequency coverage defined in + `RECEIVER_FREQUENCIES`. + + Parameters + ---------- + wavelength : `sunpy.net.attrs.Wavelength` + Input wavelength range to check + + Returns + ------- + `bool` + """ + return wavelength.min in RECEIVER_FREQUENCIES or wavelength.max in RECEIVER_FREQUENCIES + def search(self, *args, **kwargs): """ Query this client for a list of results. @@ -67,6 +79,11 @@ def search(self, *args, **kwargs): """ matchdict = self._get_match_dict(*args, **kwargs) metalist = [] + + wavelentgh = matchdict.get('Wavelength', False) + if wavelentgh and not self._check_wavelengths(wavelentgh): + return QueryResponse(metalist, client=self) + tr = TimeRange(matchdict['Start Time'], matchdict['End Time']) scraper = Scraper(self.baseurl, regex=True) filesmeta = scraper._extract_files_meta(tr, extractor=self.pattern) @@ -74,12 +91,20 @@ def search(self, *args, **kwargs): rowdict = self.post_search_hook(i, matchdict) metalist.append(rowdict) - return QueryResponse(metalist, client=self) + query_response = QueryResponse(metalist, client=self) + mask = np.full(len(query_response), True) + pol = matchdict.get('PolType') + if len(pol) == 1: + pol = pol.upper() + mask = mask & query_response['Polarisation'] == pol + + return query_response[mask] @classmethod def register_values(cls): adict = {a.Instrument: [('ILOFAR', 'Irish LOFAR STATION (IE63)')], a.Source: [('ILOFAR', 'Irish LOFAR Data Archive')], a.Provider: [('ILOFAR', 'Irish LOFAR Data Archive')], - a.Wavelength: [('*')]} + a.Wavelength: [('*')], + PolType: [('X', 'X'), ('X Linear Polarisation', 'Y Linear Polarisation')]} return adict diff --git a/radiospectra/net/sources/tests/test_ilofar.py b/radiospectra/net/sources/tests/test_ilofar.py index cc764f9..45505fb 100644 --- a/radiospectra/net/sources/tests/test_ilofar.py +++ b/radiospectra/net/sources/tests/test_ilofar.py @@ -3,15 +3,17 @@ import pytest +import astropy.units as u from sunpy.net import Fido from sunpy.net import attrs as a -from radiospectra.net.sources.ilofar import ILOFARClient +from radiospectra.net.attrs import PolType +from radiospectra.net.sources.ilofar import ILOFARMode357Client @pytest.fixture def client(): - return ILOFARClient() + return ILOFARMode357Client() @pytest.fixture @@ -43,12 +45,49 @@ def test_ilofar_client(mock_urlopen, client, html_responses): assert query[0]['Polarisation'] == 'X' +@mock.patch('sunpy.util.scraper.urlopen') +def test_ilofar_client_polarisation(mock_urlopen, client, html_responses): + mock_urlopen.return_value.read = mock.MagicMock() + mock_urlopen.return_value.read.side_effect = html_responses * 2 + mock_urlopen.close = mock.MagicMock(return_value=None) + atr = a.Time('2018/06/01', '2018/06/02') + query_x = client.search(atr, PolType('X')) + query_y = client.search(atr, PolType('Y')) + + for query, pol in zip([query_x, query_y], ['X', 'Y']): + assert len(query) == 4 + assert query[0]['Source'] == 'ILOFAR' + assert query[0]['Provider'] == 'ILOFAR' + assert query[0]['Start Time'].iso == '2018-06-01 10:00:41.000' + assert query[0]['Polarisation'] == pol + + +@mock.patch('sunpy.util.scraper.urlopen') +def test_ilofar_client_polarisation(mock_urlopen, client, html_responses): + mock_urlopen.return_value.read = mock.MagicMock() + mock_urlopen.return_value.read.side_effect = html_responses * 3 + mock_urlopen.close = mock.MagicMock(return_value=None) + atr = a.Time('2018/06/01', '2018/06/02') + query_both_low = client.search(atr, a.Wavelength(1*u.MHz, 5*u.MHz)) + query_both_high = client.search(atr, a.Wavelength(1 * u.GHz, 2 * u.GHz)) + + assert len(query_both_low) == 0 + assert len(query_both_high) == 0 + + query_low_in = client.search(atr, a.Wavelength(90*u.MHz, 1*u.GHz)) + query_both_in = client.search(atr, a.Wavelength(15*u.MHz, 230*u.MHz)) + query_high_in = client.search(atr, a.Wavelength(5*u.MHz, 90*u.MHz)) + + for query in [query_low_in, query_both_in, query_high_in]: + assert len(query) == 8 + + @pytest.mark.remote_data def test_fido(): atr = a.Time('2018/06/01', '2018/06/02') query = Fido.search(atr, a.Instrument('ILOFAR')) - assert isinstance(query[0].client, ILOFARClient) + assert isinstance(query[0].client, ILOFARMode357Client) query = query[0] assert len(query) == 8 assert query[0]['Source'] == 'ILOFAR' diff --git a/radiospectra/spectrogram2/__init__.py b/radiospectra/spectrogram2/__init__.py index f9bda63..1775518 100644 --- a/radiospectra/spectrogram2/__init__.py +++ b/radiospectra/spectrogram2/__init__.py @@ -6,7 +6,7 @@ from radiospectra.spectrogram2.sources import ( CALISTOSpectrogram, EOVSASpectrogram, - ILOFARSpectrogram, + ILOFARMode357Spectrogram, RFSSpectrogram, SWAVESSpectrogram, ) diff --git a/radiospectra/spectrogram2/sources.py b/radiospectra/spectrogram2/sources.py index 2e3ccea..42fccf3 100644 --- a/radiospectra/spectrogram2/sources.py +++ b/radiospectra/spectrogram2/sources.py @@ -5,7 +5,7 @@ from radiospectra.spectrogram2.spectrogram import GenericSpectrogram __all__ = ['SWAVESSpectrogram', 'RFSSpectrogram', 'CALISTOSpectrogram', 'EOVSASpectrogram', - 'RSTNSpectrogram', 'ILOFARSpectrogram'] + 'RSTNSpectrogram', 'ILOFARMode357Spectrogram'] class SWAVESSpectrogram(GenericSpectrogram): @@ -202,9 +202,9 @@ def is_datasource_for(cls, data, meta, **kwargs): return meta['instrument'] == 'RSTN' -class ILOFARSpectrogram(GenericSpectrogram): +class ILOFARMode357Spectrogram(GenericSpectrogram): """ - Irish LOFAR Station + Irish LOFAR Station mode 357 Spectrogram """ @property def mode(self): From f268e84d8fe2fc46ae481e7d1888853766954ffc Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Wed, 22 Sep 2021 14:26:18 +0100 Subject: [PATCH 4/8] Add addional url that is sometimes uesed --- radiospectra/net/sources/ilofar.py | 19 ++++++++++++------- radiospectra/net/sources/tests/test_ilofar.py | 18 +++++++++++++++--- 2 files changed, 27 insertions(+), 10 deletions(-) diff --git a/radiospectra/net/sources/ilofar.py b/radiospectra/net/sources/ilofar.py index 772932a..82ebd74 100644 --- a/radiospectra/net/sources/ilofar.py +++ b/radiospectra/net/sources/ilofar.py @@ -12,6 +12,8 @@ RECEIVER_FREQUENCIES = a.Wavelength(10.546875*u.MHz, 244.53125*u.MHz) +DATASET_NAMES = ['rcu357_1beam', 'rcu357_1beam_datastream'] + class ILOFARMode357Client(GenericClient): """ @@ -39,8 +41,8 @@ class ILOFARMode357Client(GenericClient): """ - baseurl = (r'https://data.lofar.ie/%Y/%m/%d/bst/kbt/rcu357_1beam/' - r'%Y%m%d_\d{6}_bst_\d{2}\S{1}.dat') + baseurl = (r'https://data.lofar.ie/%Y/%m/%d/bst/kbt/{dataset}/' + r'%Y%m%d_\d{{6}}_bst_\d{{2}}\S{{1}}.dat') pattern = r'{}/{year:4d}{month:2d}{day:2d}_{hour:2d}{minute:2d}{second:2d}' \ r'_bst_{num:2d}{Polarisation}.dat' @@ -85,11 +87,14 @@ def search(self, *args, **kwargs): return QueryResponse(metalist, client=self) tr = TimeRange(matchdict['Start Time'], matchdict['End Time']) - scraper = Scraper(self.baseurl, regex=True) - filesmeta = scraper._extract_files_meta(tr, extractor=self.pattern) - for i in filesmeta: - rowdict = self.post_search_hook(i, matchdict) - metalist.append(rowdict) + + for dataset in DATASET_NAMES: + url = self.baseurl.format(dataset=dataset) + scraper = Scraper(url, regex=True) + filesmeta = scraper._extract_files_meta(tr, extractor=self.pattern) + for i in filesmeta: + rowdict = self.post_search_hook(i, matchdict) + metalist.append(rowdict) query_response = QueryResponse(metalist, client=self) mask = np.full(len(query_response), True) diff --git a/radiospectra/net/sources/tests/test_ilofar.py b/radiospectra/net/sources/tests/test_ilofar.py index 45505fb..127edfb 100644 --- a/radiospectra/net/sources/tests/test_ilofar.py +++ b/radiospectra/net/sources/tests/test_ilofar.py @@ -30,13 +30,15 @@ def html_responses(): @mock.patch('sunpy.util.scraper.urlopen') def test_ilofar_client(mock_urlopen, client, html_responses): mock_urlopen.return_value.read = mock.MagicMock() - mock_urlopen.return_value.read.side_effect = html_responses + mock_urlopen.return_value.read.side_effect = html_responses * 2 mock_urlopen.close = mock.MagicMock(return_value=None) atr = a.Time('2018/06/01', '2018/06/02') query = client.search(atr) called_urls = ['https://data.lofar.ie/2018/06/01/bst/kbt/rcu357_1beam/', - 'https://data.lofar.ie/2018/06/02/bst/kbt/rcu357_1beam/'] + 'https://data.lofar.ie/2018/06/02/bst/kbt/rcu357_1beam/', + 'https://data.lofar.ie/2018/06/01/bst/kbt/rcu357_1beam_datastream/', + 'https://data.lofar.ie/2018/06/02/bst/kbt/rcu357_1beam_datastream/'] assert called_urls == [call[0][0] for call in mock_urlopen.call_args_list] assert len(query) == 8 assert query[0]['Source'] == 'ILOFAR' @@ -65,7 +67,7 @@ def test_ilofar_client_polarisation(mock_urlopen, client, html_responses): @mock.patch('sunpy.util.scraper.urlopen') def test_ilofar_client_polarisation(mock_urlopen, client, html_responses): mock_urlopen.return_value.read = mock.MagicMock() - mock_urlopen.return_value.read.side_effect = html_responses * 3 + mock_urlopen.return_value.read.side_effect = html_responses * 6 mock_urlopen.close = mock.MagicMock(return_value=None) atr = a.Time('2018/06/01', '2018/06/02') query_both_low = client.search(atr, a.Wavelength(1*u.MHz, 5*u.MHz)) @@ -94,3 +96,13 @@ def test_fido(): assert query[0]['Provider'] == 'ILOFAR' assert query[0]['Start Time'].iso == '2018-06-01 10:00:41.000' assert query[0]['Polarisation'] == 'X' + + +@pytest.mark.remote_data +def test_fido_other_dataset(): + atr = a.Time('2021/08/01', '2021/10/01') + query = Fido.search(atr, a.Instrument('ILOFAR')) + + assert isinstance(query[0].client, ILOFARMode357Client) + query = query[0] + assert len(query) == 34 From 66db8c13021f599da9354c95b12311d6b884df93 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Wed, 22 Sep 2021 22:56:45 +0100 Subject: [PATCH 5/8] Fix bug introduced by update for ILOFAR --- radiospectra/net/sources/ilofar.py | 26 +++++++++++++++--------- radiospectra/spectrogram2/sources.py | 2 +- radiospectra/spectrogram2/spectrogram.py | 5 ++++- 3 files changed, 21 insertions(+), 12 deletions(-) diff --git a/radiospectra/net/sources/ilofar.py b/radiospectra/net/sources/ilofar.py index 82ebd74..36a0f4d 100644 --- a/radiospectra/net/sources/ilofar.py +++ b/radiospectra/net/sources/ilofar.py @@ -24,20 +24,26 @@ class ILOFARMode357Client(GenericClient): -------- >>> import radiospectra.net >>> from sunpy.net import Fido, attrs as a - >>> results = Fido.search(a.Time("2010/10/01", "2010/10/02"), + >>> results = Fido.search(a.Time("2021/09/01", "2021/09/21"), ... a.Instrument('ILOFAR')) # doctest: +REMOTE_DATA >>> results #doctest: +REMOTE_DATA - + - 4 Results from the WAVESClient: - Start Time End Time ... Provider Wavelength [2] - ... kHz - ----------------------- ----------------------- ... -------- ----------------- - 2010-10-01 00:00:00.000 2010-10-01 23:59:59.999 ... NASA 20.0 .. 1040.0 - 2010-10-02 00:00:00.000 2010-10-02 23:59:59.999 ... NASA 20.0 .. 1040.0 - 2010-10-01 00:00:00.000 2010-10-01 23:59:59.999 ... NASA 1075.0 .. 13825.0 - 2010-10-02 00:00:00.000 2010-10-02 23:59:59.999 ... NASA 1075.0 .. 13825.0 + 10 Results from the ILOFARMode357Client: + Start Time End Time Instrument ... num Polarisation + ----------------------- ----------------------- ---------- ... --- ------------ + 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ILOFAR ... 0 X + 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ILOFAR ... 0 Y + 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ILOFAR ... 0 X + 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ILOFAR ... 0 Y + 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ILOFAR ... 0 X + 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ILOFAR ... 0 Y + 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ILOFAR ... 0 X + 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ILOFAR ... 0 Y + 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ILOFAR ... 0 X + 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ILOFAR ... 0 Y + """ diff --git a/radiospectra/spectrogram2/sources.py b/radiospectra/spectrogram2/sources.py index 42fccf3..69a95b9 100644 --- a/radiospectra/spectrogram2/sources.py +++ b/radiospectra/spectrogram2/sources.py @@ -5,7 +5,7 @@ from radiospectra.spectrogram2.spectrogram import GenericSpectrogram __all__ = ['SWAVESSpectrogram', 'RFSSpectrogram', 'CALISTOSpectrogram', 'EOVSASpectrogram', - 'RSTNSpectrogram', 'ILOFARMode357Spectrogram'] + 'RSTNSpectrogram', 'ILOFARMode357Spectrogram', 'WAVESSpectrogram'] class SWAVESSpectrogram(GenericSpectrogram): diff --git a/radiospectra/spectrogram2/spectrogram.py b/radiospectra/spectrogram2/spectrogram.py index 1ead98a..8e03d68 100644 --- a/radiospectra/spectrogram2/spectrogram.py +++ b/radiospectra/spectrogram2/spectrogram.py @@ -115,7 +115,10 @@ def parse_path(path, f, **kwargs): raise ValueError('path must be a pathlib.Path object') path = path.expanduser() if is_file(path): - return list(f(path, **kwargs)) + read_file = f(path, **kwargs) + if isinstance(read_file, tuple): + read_file = [read_file] + return read_file elif is_dir(path): read_files = [] for afile in sorted(path.glob('*')): From a919b78fa55d24cc1e3bbf960087ad50bdafd343 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Tue, 12 Oct 2021 13:23:45 +0200 Subject: [PATCH 6/8] Fix issue with LOFAR client and spectrogram * BST files with 01 at the end are not valid so ignore in fido client * Remove duble entry of polarisation data from client results * BST files can be corrupt, try to detect and use intact part of file in spectrogram --- radiospectra/net/sources/ilofar.py | 32 ++++++++++--------- radiospectra/net/sources/tests/test_ilofar.py | 2 +- radiospectra/spectrogram2/spectrogram.py | 6 ++++ 3 files changed, 24 insertions(+), 16 deletions(-) diff --git a/radiospectra/net/sources/ilofar.py b/radiospectra/net/sources/ilofar.py index 36a0f4d..729f017 100644 --- a/radiospectra/net/sources/ilofar.py +++ b/radiospectra/net/sources/ilofar.py @@ -30,28 +30,27 @@ class ILOFARMode357Client(GenericClient): - 10 Results from the ILOFARMode357Client: - Start Time End Time Instrument ... num Polarisation - ----------------------- ----------------------- ---------- ... --- ------------ - 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ILOFAR ... 0 X - 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ILOFAR ... 0 Y - 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ILOFAR ... 0 X - 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ILOFAR ... 0 Y - 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ILOFAR ... 0 X - 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ILOFAR ... 0 Y - 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ILOFAR ... 0 X - 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ILOFAR ... 0 Y - 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ILOFAR ... 0 X - 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ILOFAR ... 0 Y + Start Time End Time ... PolType Polarisation + ----------------------- ----------------------- ... ------- ------------ + 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ... X X + 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ... X Y + 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ... X X + 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ... X Y + 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ... X X + 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ... X Y + 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ... X X + 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ... X Y + 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ... X X + 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ... X Y """ baseurl = (r'https://data.lofar.ie/%Y/%m/%d/bst/kbt/{dataset}/' - r'%Y%m%d_\d{{6}}_bst_\d{{2}}\S{{1}}.dat') + r'%Y%m%d_\d{{6}}_bst_00\S{{1}}.dat') pattern = r'{}/{year:4d}{month:2d}{day:2d}_{hour:2d}{minute:2d}{second:2d}' \ - r'_bst_{num:2d}{Polarisation}.dat' + r'_bst_00{Polarisation}.dat' @classmethod def _check_wavelengths(cls, wavelength): @@ -109,6 +108,9 @@ def search(self, *args, **kwargs): pol = pol.upper() mask = mask & query_response['Polarisation'] == pol + if query_response: + query_response.remove_column('PolType') + return query_response[mask] @classmethod diff --git a/radiospectra/net/sources/tests/test_ilofar.py b/radiospectra/net/sources/tests/test_ilofar.py index 127edfb..4e7c677 100644 --- a/radiospectra/net/sources/tests/test_ilofar.py +++ b/radiospectra/net/sources/tests/test_ilofar.py @@ -105,4 +105,4 @@ def test_fido_other_dataset(): assert isinstance(query[0].client, ILOFARMode357Client) query = query[0] - assert len(query) == 34 + assert len(query) == 38 diff --git a/radiospectra/spectrogram2/spectrogram.py b/radiospectra/spectrogram2/spectrogram.py index 8e03d68..5dea66a 100644 --- a/radiospectra/spectrogram2/spectrogram.py +++ b/radiospectra/spectrogram2/spectrogram.py @@ -20,6 +20,7 @@ from astropy.io.fits import Header from astropy.time import Time from astropy.visualization import quantity_support +from sunpy import log from sunpy.data import cache from sunpy.io import fits from sunpy.net import attrs as a @@ -567,6 +568,11 @@ def subband_to_freq(subband, obs_mode): polarisation = file.stem[-1] num_times = data.shape[0] / num_subbands + if not num_times.is_integer(): + log.warning('BST file seems incomplete dropping incomplete frequencies') + num_times = np.floor(num_times).astype(int) + truncate = num_times * num_subbands + data = data[:truncate] data = data.reshape(-1, num_subbands).T # (Freq x Time).T = (Time x Freq) dt = np.arange(num_times) * 1 * u.s start_time = Time.strptime(file.name.split('_bst')[0], "%Y%m%d_%H%M%S") From 7da679732315bf233ff3dfdc507fbf9936502227 Mon Sep 17 00:00:00 2001 From: Nabil Freij Date: Thu, 10 Nov 2022 18:40:51 -0800 Subject: [PATCH 7/8] fix bad merge --- CHANGELOG.rst | 4 + changelog/57.feature.rst | 1 - radiospectra/mixins.py | 5 +- radiospectra/net/__init__.py | 11 +- radiospectra/net/sources/ilofar.py | 35 +- radiospectra/net/sources/tests/test_ilofar.py | 72 +- radiospectra/spectrogram/sources/ilofar357.py | 23 + .../sources/tests/test_ilofar357.py | 28 + .../spectrogram/spectrogram_factory.py | 46 ++ radiospectra/spectrogram2/__init__.py | 13 - radiospectra/spectrogram2/sources.py | 219 ----- radiospectra/spectrogram2/spectrogram.py | 773 ------------------ .../spectrogram2/tests/test_ilofar.py | 28 - radiospectra/utils.py | 29 + 14 files changed, 197 insertions(+), 1090 deletions(-) delete mode 100644 changelog/57.feature.rst create mode 100644 radiospectra/spectrogram/sources/ilofar357.py create mode 100644 radiospectra/spectrogram/sources/tests/test_ilofar357.py delete mode 100644 radiospectra/spectrogram2/__init__.py delete mode 100644 radiospectra/spectrogram2/sources.py delete mode 100644 radiospectra/spectrogram2/spectrogram.py delete mode 100644 radiospectra/spectrogram2/tests/test_ilofar.py create mode 100644 radiospectra/utils.py diff --git a/CHANGELOG.rst b/CHANGELOG.rst index de8ffc4..b64d1f8 100644 --- a/CHANGELOG.rst +++ b/CHANGELOG.rst @@ -8,6 +8,10 @@ Breaking Changes - The new ``Spectrogram2`` class has been renamed to ``Spectrogram``. (`#76 `__) - Adding colorbar functionality to ``plot`` (`#80 `__) +Features +-------- +- Add `sunpy.net.Fido` client `~radiospectra.net.sources.ilofar.ILOFARMode357` and spectrogram class `~radiospectra.spectrogram2.sources.ILOFARMode357` for ILOFAR mode 357 observations. (`#57 `__) + 0.4.0 (2022-05-24) ================== diff --git a/changelog/57.feature.rst b/changelog/57.feature.rst deleted file mode 100644 index 58b9927..0000000 --- a/changelog/57.feature.rst +++ /dev/null @@ -1 +0,0 @@ -Add `sunpy.net.Fido` client `~radiospectra.net.sources.ilofar.ILOFARMode357` and spectrogram class `~radiospectra.spectrogram2.sources.ILOFARMode357` for ILOFAR mode 357 observations. diff --git a/radiospectra/mixins.py b/radiospectra/mixins.py index 05d44b0..6ecdfc5 100644 --- a/radiospectra/mixins.py +++ b/radiospectra/mixins.py @@ -37,7 +37,10 @@ def plot(self, axes=None, **kwargs): axes.set_title(title) axes.plot(self.times.datetime[[0, -1]], self.frequencies[[0, -1]], linestyle="None", marker="None") - ret = axes.pcolormesh(self.times.datetime, self.frequencies.value, data[:-1, :-1], shading="auto", **kwargs) + if self.times.shape[0] == self.data.shape[0] and self.frequencies.shape[0] == self.data.shape[1]: + ret = axes.pcolormesh(self.times.datetime, self.frequencies.value, data, shading="auto", **kwargs) + else: + ret = axes.pcolormesh(self.times.datetime, self.frequencies.value, data[:-1, :-1], shading="auto", **kwargs) axes.set_xlim(self.times.datetime[0], self.times.datetime[-1]) locator = mdates.AutoDateLocator(minticks=4, maxticks=8) formatter = mdates.ConciseDateFormatter(locator) diff --git a/radiospectra/net/__init__.py b/radiospectra/net/__init__.py index 1206211..886ae3d 100644 --- a/radiospectra/net/__init__.py +++ b/radiospectra/net/__init__.py @@ -7,5 +7,12 @@ from radiospectra.net.sources.stereo import SWAVESClient from radiospectra.net.sources.wind import WAVESClient -__all__ = ['CALLISTOClient', 'EOVSAClient', 'RFSClient', 'SWAVESClient', 'RSTNClient', - 'WAVESClient', 'ILOFARMode357Client'] +__all__ = [ + "CALLISTOClient", + "EOVSAClient", + "RFSClient", + "SWAVESClient", + "RSTNClient", + "WAVESClient", + "ILOFARMode357Client", +] diff --git a/radiospectra/net/sources/ilofar.py b/radiospectra/net/sources/ilofar.py index 729f017..633cf54 100644 --- a/radiospectra/net/sources/ilofar.py +++ b/radiospectra/net/sources/ilofar.py @@ -6,13 +6,12 @@ from sunpy.time import TimeRange from sunpy.util.scraper import Scraper -__all__ = ['ILOFARMode357Client'] - from radiospectra.net.attrs import PolType -RECEIVER_FREQUENCIES = a.Wavelength(10.546875*u.MHz, 244.53125*u.MHz) +__all__ = ["ILOFARMode357Client"] -DATASET_NAMES = ['rcu357_1beam', 'rcu357_1beam_datastream'] +RECEIVER_FREQUENCIES = a.Wavelength(10.546875 * u.MHz, 244.53125 * u.MHz) +DATASET_NAMES = ["rcu357_1beam", "rcu357_1beam_datastream"] class ILOFARMode357Client(GenericClient): @@ -46,11 +45,9 @@ class ILOFARMode357Client(GenericClient): """ - baseurl = (r'https://data.lofar.ie/%Y/%m/%d/bst/kbt/{dataset}/' - r'%Y%m%d_\d{{6}}_bst_00\S{{1}}.dat') + baseurl = r"https://data.lofar.ie/%Y/%m/%d/bst/kbt/{dataset}/" r"%Y%m%d_\d{{6}}_bst_00\S{{1}}.dat" - pattern = r'{}/{year:4d}{month:2d}{day:2d}_{hour:2d}{minute:2d}{second:2d}' \ - r'_bst_00{Polarisation}.dat' + pattern = r"{}/{year:4d}{month:2d}{day:2d}_{hour:2d}{minute:2d}{second:2d}" r"_bst_00{Polarisation}.dat" @classmethod def _check_wavelengths(cls, wavelength): @@ -87,11 +84,11 @@ def search(self, *args, **kwargs): matchdict = self._get_match_dict(*args, **kwargs) metalist = [] - wavelentgh = matchdict.get('Wavelength', False) + wavelentgh = matchdict.get("Wavelength", False) if wavelentgh and not self._check_wavelengths(wavelentgh): return QueryResponse(metalist, client=self) - tr = TimeRange(matchdict['Start Time'], matchdict['End Time']) + tr = TimeRange(matchdict["Start Time"], matchdict["End Time"]) for dataset in DATASET_NAMES: url = self.baseurl.format(dataset=dataset) @@ -103,21 +100,23 @@ def search(self, *args, **kwargs): query_response = QueryResponse(metalist, client=self) mask = np.full(len(query_response), True) - pol = matchdict.get('PolType') + pol = matchdict.get("PolType") if len(pol) == 1: pol = pol.upper() - mask = mask & query_response['Polarisation'] == pol + mask = mask & query_response["Polarisation"] == pol if query_response: - query_response.remove_column('PolType') + query_response.remove_column("PolType") return query_response[mask] @classmethod def register_values(cls): - adict = {a.Instrument: [('ILOFAR', 'Irish LOFAR STATION (IE63)')], - a.Source: [('ILOFAR', 'Irish LOFAR Data Archive')], - a.Provider: [('ILOFAR', 'Irish LOFAR Data Archive')], - a.Wavelength: [('*')], - PolType: [('X', 'X'), ('X Linear Polarisation', 'Y Linear Polarisation')]} + adict = { + a.Instrument: [("ILOFAR", "Irish LOFAR STATION (IE63)")], + a.Source: [("ILOFAR", "Irish LOFAR Data Archive")], + a.Provider: [("ILOFAR", "Irish LOFAR Data Archive")], + a.Wavelength: [("*")], + PolType: [("X", "X"), ("X Linear Polarisation", "Y Linear Polarisation")], + } return adict diff --git a/radiospectra/net/sources/tests/test_ilofar.py b/radiospectra/net/sources/tests/test_ilofar.py index 4e7c677..6535cda 100644 --- a/radiospectra/net/sources/tests/test_ilofar.py +++ b/radiospectra/net/sources/tests/test_ilofar.py @@ -18,67 +18,69 @@ def client(): @pytest.fixture def html_responses(): - paths = [Path(__file__).parent / 'data' / n for n in ['ilofar_resp1.html', 'ilofar_resp2.html']] + paths = [Path(__file__).parent / "data" / n for n in ["ilofar_resp1.html", "ilofar_resp2.html"]] response_htmls = [] for p in paths: - with p.open('r') as f: + with p.open("r") as f: response_htmls.append(f.read()) return response_htmls -@mock.patch('sunpy.util.scraper.urlopen') +@mock.patch("sunpy.util.scraper.urlopen") def test_ilofar_client(mock_urlopen, client, html_responses): mock_urlopen.return_value.read = mock.MagicMock() mock_urlopen.return_value.read.side_effect = html_responses * 2 mock_urlopen.close = mock.MagicMock(return_value=None) - atr = a.Time('2018/06/01', '2018/06/02') + atr = a.Time("2018/06/01", "2018/06/02") query = client.search(atr) - called_urls = ['https://data.lofar.ie/2018/06/01/bst/kbt/rcu357_1beam/', - 'https://data.lofar.ie/2018/06/02/bst/kbt/rcu357_1beam/', - 'https://data.lofar.ie/2018/06/01/bst/kbt/rcu357_1beam_datastream/', - 'https://data.lofar.ie/2018/06/02/bst/kbt/rcu357_1beam_datastream/'] + called_urls = [ + "https://data.lofar.ie/2018/06/01/bst/kbt/rcu357_1beam/", + "https://data.lofar.ie/2018/06/02/bst/kbt/rcu357_1beam/", + "https://data.lofar.ie/2018/06/01/bst/kbt/rcu357_1beam_datastream/", + "https://data.lofar.ie/2018/06/02/bst/kbt/rcu357_1beam_datastream/", + ] assert called_urls == [call[0][0] for call in mock_urlopen.call_args_list] assert len(query) == 8 - assert query[0]['Source'] == 'ILOFAR' - assert query[0]['Provider'] == 'ILOFAR' - assert query[0]['Start Time'].iso == '2018-06-01 10:00:41.000' - assert query[0]['Polarisation'] == 'X' + assert query[0]["Source"] == "ILOFAR" + assert query[0]["Provider"] == "ILOFAR" + assert query[0]["Start Time"].iso == "2018-06-01 10:00:41.000" + assert query[0]["Polarisation"] == "X" -@mock.patch('sunpy.util.scraper.urlopen') +@mock.patch("sunpy.util.scraper.urlopen") def test_ilofar_client_polarisation(mock_urlopen, client, html_responses): mock_urlopen.return_value.read = mock.MagicMock() mock_urlopen.return_value.read.side_effect = html_responses * 2 mock_urlopen.close = mock.MagicMock(return_value=None) - atr = a.Time('2018/06/01', '2018/06/02') - query_x = client.search(atr, PolType('X')) - query_y = client.search(atr, PolType('Y')) + atr = a.Time("2018/06/01", "2018/06/02") + query_x = client.search(atr, PolType("X")) + query_y = client.search(atr, PolType("Y")) - for query, pol in zip([query_x, query_y], ['X', 'Y']): + for query, pol in zip([query_x, query_y], ["X", "Y"]): assert len(query) == 4 - assert query[0]['Source'] == 'ILOFAR' - assert query[0]['Provider'] == 'ILOFAR' - assert query[0]['Start Time'].iso == '2018-06-01 10:00:41.000' - assert query[0]['Polarisation'] == pol + assert query[0]["Source"] == "ILOFAR" + assert query[0]["Provider"] == "ILOFAR" + assert query[0]["Start Time"].iso == "2018-06-01 10:00:41.000" + assert query[0]["Polarisation"] == pol -@mock.patch('sunpy.util.scraper.urlopen') +@mock.patch("sunpy.util.scraper.urlopen") def test_ilofar_client_polarisation(mock_urlopen, client, html_responses): mock_urlopen.return_value.read = mock.MagicMock() mock_urlopen.return_value.read.side_effect = html_responses * 6 mock_urlopen.close = mock.MagicMock(return_value=None) - atr = a.Time('2018/06/01', '2018/06/02') - query_both_low = client.search(atr, a.Wavelength(1*u.MHz, 5*u.MHz)) + atr = a.Time("2018/06/01", "2018/06/02") + query_both_low = client.search(atr, a.Wavelength(1 * u.MHz, 5 * u.MHz)) query_both_high = client.search(atr, a.Wavelength(1 * u.GHz, 2 * u.GHz)) assert len(query_both_low) == 0 assert len(query_both_high) == 0 - query_low_in = client.search(atr, a.Wavelength(90*u.MHz, 1*u.GHz)) - query_both_in = client.search(atr, a.Wavelength(15*u.MHz, 230*u.MHz)) - query_high_in = client.search(atr, a.Wavelength(5*u.MHz, 90*u.MHz)) + query_low_in = client.search(atr, a.Wavelength(90 * u.MHz, 1 * u.GHz)) + query_both_in = client.search(atr, a.Wavelength(15 * u.MHz, 230 * u.MHz)) + query_high_in = client.search(atr, a.Wavelength(5 * u.MHz, 90 * u.MHz)) for query in [query_low_in, query_both_in, query_high_in]: assert len(query) == 8 @@ -86,22 +88,22 @@ def test_ilofar_client_polarisation(mock_urlopen, client, html_responses): @pytest.mark.remote_data def test_fido(): - atr = a.Time('2018/06/01', '2018/06/02') - query = Fido.search(atr, a.Instrument('ILOFAR')) + atr = a.Time("2018/06/01", "2018/06/02") + query = Fido.search(atr, a.Instrument("ILOFAR")) assert isinstance(query[0].client, ILOFARMode357Client) query = query[0] assert len(query) == 8 - assert query[0]['Source'] == 'ILOFAR' - assert query[0]['Provider'] == 'ILOFAR' - assert query[0]['Start Time'].iso == '2018-06-01 10:00:41.000' - assert query[0]['Polarisation'] == 'X' + assert query[0]["Source"] == "ILOFAR" + assert query[0]["Provider"] == "ILOFAR" + assert query[0]["Start Time"].iso == "2018-06-01 10:00:41.000" + assert query[0]["Polarisation"] == "X" @pytest.mark.remote_data def test_fido_other_dataset(): - atr = a.Time('2021/08/01', '2021/10/01') - query = Fido.search(atr, a.Instrument('ILOFAR')) + atr = a.Time("2021/08/01", "2021/10/01") + query = Fido.search(atr, a.Instrument("ILOFAR")) assert isinstance(query[0].client, ILOFARMode357Client) query = query[0] diff --git a/radiospectra/spectrogram/sources/ilofar357.py b/radiospectra/spectrogram/sources/ilofar357.py new file mode 100644 index 0000000..1dd2273 --- /dev/null +++ b/radiospectra/spectrogram/sources/ilofar357.py @@ -0,0 +1,23 @@ +from radiospectra.spectrogram.spectrogrambase import GenericSpectrogram + +__all__ = [ + "ILOFARMode357Spectrogram", +] + + +class ILOFARMode357Spectrogram(GenericSpectrogram): + """ + Irish LOFAR Station mode 357 Spectrogram + """ + + @property + def mode(self): + return self.meta.get("mode") + + @property + def polarisation(self): + return self.meta.get("polarisation") + + @classmethod + def is_datasource_for(cls, data, meta, **kwargs): + return meta["instrument"] == "ILOFAR" diff --git a/radiospectra/spectrogram/sources/tests/test_ilofar357.py b/radiospectra/spectrogram/sources/tests/test_ilofar357.py new file mode 100644 index 0000000..9d42b60 --- /dev/null +++ b/radiospectra/spectrogram/sources/tests/test_ilofar357.py @@ -0,0 +1,28 @@ +from pathlib import Path +from unittest import mock + +import numpy as np + +from radiospectra.spectrogram import Spectrogram + + +@mock.patch("numpy.fromfile") +@mock.patch("radiospectra.spectrogram2.spectrogram.is_file") +def test_ilofar(mock_is_file, mock_fromfile): + mock_fromfile.return_value = np.ones(10117216) + mock_is_file.return_value = True + + spec = Spectrogram(Path("20180602_063247_bst_00X.dat")) + assert len(spec) == 3 + assert spec[0].polarisation == "X" + assert spec[0].start_time.iso == "2018-06-02 06:32:47.000" + assert spec[0].end_time.iso == "2018-06-02 12:18:18.000" + assert spec[0].mode == 3 + assert spec[0].frequencies[0].to_value("MHz") == 10.546875 + assert spec[0].frequencies[-1].to_value("MHz") == 88.28125 + assert spec[1].mode == 5 + assert spec[1].frequencies[0].to_value("MHz") == 110.546875 + assert spec[1].frequencies[-1].to_value("MHz") == 188.28125 + assert spec[2].mode == 7 + assert spec[2].frequencies[0].to_value("MHz") == 210.546875 + assert spec[2].frequencies[-1].to_value("MHz") == 244.53125 diff --git a/radiospectra/spectrogram/spectrogram_factory.py b/radiospectra/spectrogram/spectrogram_factory.py index 80f5554..efb7d06 100644 --- a/radiospectra/spectrogram/spectrogram_factory.py +++ b/radiospectra/spectrogram/spectrogram_factory.py @@ -16,6 +16,7 @@ from astropy.io import fits from astropy.io.fits import Header from astropy.time import Time +from sunpy import log from sunpy.data import cache from sunpy.net import attrs as a from sunpy.time import parse_time @@ -32,6 +33,7 @@ from radiospectra.exceptions import NoSpectrogramInFileError, SpectraMetaValidationError from radiospectra.spectrogram.spectrogrambase import GenericSpectrogram +from radiospectra.utils import subband_to_freq SUPPORTED_ARRAY_TYPES = (np.ndarray,) try: @@ -273,6 +275,50 @@ def _read_dat(file): meta["times"] = meta["start_time"] + times meta["end_time"] = meta["start_time"] + times[-1] return data, meta + elif "bst" in file.name: + subbands = (np.arange(54, 454, 2), np.arange(54, 454, 2), np.arange(54, 230, 2)) + num_subbands = 488 + + data = np.fromfile(file) + polarisation = file.stem[-1] + + num_times = data.shape[0] / num_subbands + if not num_times.is_integer(): + log.warning("BST file seems incomplete dropping incomplete frequencies") + num_times = np.floor(num_times).astype(int) + truncate = num_times * num_subbands + data = data[:truncate] + data = data.reshape(-1, num_subbands).T # (Freq x Time).T = (Time x Freq) + dt = np.arange(num_times) * 1 * u.s + start_time = Time.strptime(file.name.split("_bst")[0], "%Y%m%d_%H%M%S") + times = start_time + dt + + obs_mode = (3, 5, 7) + + freqs = [subband_to_freq(sb, mode) for sb, mode in zip(subbands, obs_mode)] + + # 1st 200 sbs mode 3, next 200 sbs mode 5, last 88 sbs mode 7 + spec = {0: data[:200, :], 1: data[200:400, :], 2: data[400:, :]} + data_header_pairs = [] + for i in range(3): + meta = { + "instrument": "ILOFAR", + "observatory": "Birr (IE613)", + "start_time": times[0], + "mode": obs_mode[i], + "wavelength": a.Wavelength(freqs[i][0], freqs[i][-1]), + "freqs": freqs[i], + "times": times, + "end_time": times[-1], + "detector": "ILOFAR", + "polarisation": polarisation, + } + + data_header_pairs.append((spec[i], meta)) + + return data_header_pairs + else: + raise ValueError(f"File {file} not supported.") @staticmethod def _read_srs(file): diff --git a/radiospectra/spectrogram2/__init__.py b/radiospectra/spectrogram2/__init__.py deleted file mode 100644 index 1775518..0000000 --- a/radiospectra/spectrogram2/__init__.py +++ /dev/null @@ -1,13 +0,0 @@ -""" -The spectrogram2 module aims to provide a more sunpy-like interface to radiospectra data similar to -that of `~sunpy.map.Map` and `~sunpy.timeseries.TimeSeries`. -""" - -from radiospectra.spectrogram2.sources import ( - CALISTOSpectrogram, - EOVSASpectrogram, - ILOFARMode357Spectrogram, - RFSSpectrogram, - SWAVESSpectrogram, -) -from radiospectra.spectrogram2.spectrogram import Spectrogram diff --git a/radiospectra/spectrogram2/sources.py b/radiospectra/spectrogram2/sources.py deleted file mode 100644 index 69a95b9..0000000 --- a/radiospectra/spectrogram2/sources.py +++ /dev/null @@ -1,219 +0,0 @@ - -import astropy.units as u -from astropy.coordinates.earth import EarthLocation - -from radiospectra.spectrogram2.spectrogram import GenericSpectrogram - -__all__ = ['SWAVESSpectrogram', 'RFSSpectrogram', 'CALISTOSpectrogram', 'EOVSASpectrogram', - 'RSTNSpectrogram', 'ILOFARMode357Spectrogram', 'WAVESSpectrogram'] - - -class SWAVESSpectrogram(GenericSpectrogram): - """ - STEREO Waves or S/WAVES, SWAVES Spectrogram - - Examples - -------- - >>> import radiospectra.net - >>> from sunpy.net import Fido, attrs as a - >>> from radiospectra.spectrogram2 import Spectrogram - >>> from radiospectra.net import attrs as ra - >>> query = Fido.search(a.Time('2019/10/05 23:00', '2019/10/06 00:59'), #doctest: +REMOTE_DATA - ... a.Instrument('SWAVES')) #doctest: +REMOTE_DATA - >>> downloaded = Fido.fetch(query[1][0]) #doctest: +REMOTE_DATA - >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA - >>> spec #doctest: +REMOTE_DATA - - >>> spec.plot() #doctest: +REMOTE_DATA - """ - def __init__(self, data, meta, **kwargs): - super().__init__(meta=meta, data=data, **kwargs) - - @property - def receiver(self): - """ - The name of the receiver - """ - return self.meta['receiver'] - - @classmethod - def is_datasource_for(cls, data, meta, **kwargs): - return meta['instrument'] == 'swaves' - - -class WAVESSpectrogram(GenericSpectrogram): - """ - Wind Waves Spectrogram - - Examples - -------- - >>> import radiospectra.net - >>> from sunpy.net import Fido, attrs as a - >>> from radiospectra.spectrogram2 import Spectrogram - >>> from radiospectra.net import attrs as ra - >>> query = Fido.search(a.Time('2019/10/05 23:00', '2019/10/06 00:59'), #doctest: +REMOTE_DATA - ... a.Instrument('WAVES')) #doctest: +REMOTE_DATA - >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA - >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA - >>> spec #doctest: +REMOTE_DATA - - >>> spec.plot() #doctest: +REMOTE_DATA - """ - def __init__(self, data, meta, **kwargs): - super().__init__(meta=meta, data=data, **kwargs) - - @property - def receiver(self): - """ - The name of the receiver - """ - return self.meta['receiver'] - - @property - def background(self): - """ - The background subtracted from the data - """ - return self.meta.bg - - @classmethod - def is_datasource_for(cls, data, meta, **kwargs): - return meta.get('instrument', None) == 'WAVES' - - -class RFSSpectrogram(GenericSpectrogram): - """ - Parker Solar Probe FIELDS/Radio Frequency Spectrometer (RFS) Spectrogram - - >>> import radiospectra.net - >>> from sunpy.net import Fido, attrs as a - >>> from radiospectra.spectrogram2 import Spectrogram - >>> from radiospectra.net import attrs as ra - >>> query = Fido.search(a.Time('2019/10/05 23:00', '2019/10/06 00:59'), #doctest: +REMOTE_DATA - ... a.Instrument.rfs) #doctest: +REMOTE_DATA - >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA - >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA - >>> spec #doctest: +REMOTE_DATA - - >>> spec.plot() #doctest: +REMOTE_DATA - """ - def __init__(self, data, meta, **kwargs): - super().__init__(meta=meta, data=data, **kwargs) - - @property - def level(self): - return self.meta['cdf_meta']['Data_type'].split('>')[0] - - @property - def version(self): - return int(self.meta['cdf_meta']['Data_version']) - - @classmethod - def is_datasource_for(cls, data, meta, **kwargs): - return (meta['observatory'] == 'PSP' and meta['instrument'] == 'FIELDS/RFS' - and meta['detector'] in ('lfr', 'hfr')) - - -class CALISTOSpectrogram(GenericSpectrogram): - """ - CALISTO Spectrogram from the e-CALISTO network - - Examples - -------- - >>> import radiospectra.net - >>> from sunpy.net import Fido, attrs as a - >>> from radiospectra.spectrogram2 import Spectrogram - >>> from radiospectra.net import attrs as ra - >>> query = Fido.search(a.Time('2019/10/05 23:00', '2019/10/06 00:59'), #doctest: +REMOTE_DATA - ... a.Instrument('eCALLISTO'), ra.Observatory('ALASKA')) #doctest: +REMOTE_DATA - >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA - >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA - >>> spec #doctest: +REMOTE_DATA - - >>> spec.plot() #doctest: +REMOTE_DATA - """ - def __init__(self, data, meta, **kwargs): - super().__init__(meta=meta, data=data, **kwargs) - - @property - def observatory_location(self): - lat = self.meta['fits_meta']['OBS_LAT'] * u.deg * 1.0 if \ - self.meta['fits_meta']['OBS_LAC'] == 'N' else -1.0 - lon = self.meta['fits_meta']['OBS_LON'] * u.deg * 1.0 if \ - self.meta['fits_meta']['OBS_LOC'] == 'E' else -1.0 - height = self.meta['fits_meta']['OBS_ALT'] * u.m - return EarthLocation(lat=lat, lon=lon, height=height) - - @classmethod - def is_datasource_for(cls, data, meta, **kwargs): - return meta['instrument'] == 'e-CALLISTO' or meta['detector'] == 'e-CALLISTO' - - -class EOVSASpectrogram(GenericSpectrogram): - """ - Extend Owen Valley Array (EOVSA) Spectrogram - - Examples - -------- - >>> import radiospectra.net - >>> from sunpy.net import Fido, attrs as a - >>> from radiospectra.spectrogram2 import Spectrogram - >>> query = Fido.search(a.Time('2021/05/07 00:00', '2021/05/07 23:00'), a.Instrument.eovsa) #doctest: +REMOTE_DATA - >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA - >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA - >>> spec #doctest: +REMOTE_DATA - - >>> spec.plot() #doctest: +REMOTE_DATA - """ - def __init__(self, data, meta, **kwargs): - super().__init__(meta=meta, data=data, **kwargs) - - @property - def polarisation(self): - return self.meta['fits_meta']['POLARIZA'] - - @classmethod - def is_datasource_for(cls, data, meta, **kwargs): - return meta['instrument'] == 'EOVSA' or meta['detector'] == 'EOVSA' - - # TODO fix time gaps for plots need to render them as gaps - # can prob do when generateing proper pcolormesh gird but then prob doesn't belong here - - -class RSTNSpectrogram(GenericSpectrogram): - """ - Radio Solar Telescope Network - - Examples - -------- - >>> import radiospectra.net - >>> from sunpy.net import Fido, attrs as a - >>> from radiospectra.spectrogram2 import Spectrogram - >>> query = Fido.search(a.Time('2017/09/07 00:00', '2017/09/07 23:00'), a.Instrument.rstn) #doctest: +REMOTE_DATA - >>> downloaded = Fido.fetch(query[0][0]) #doctest: +REMOTE_DATA - >>> spec = Spectrogram(downloaded[0]) #doctest: +REMOTE_DATA - >>> spec #doctest: +REMOTE_DATA - - >>> spec.plot() #doctest: +REMOTE_DATA - """ - - @classmethod - def is_datasource_for(cls, data, meta, **kwargs): - return meta['instrument'] == 'RSTN' - - -class ILOFARMode357Spectrogram(GenericSpectrogram): - """ - Irish LOFAR Station mode 357 Spectrogram - """ - @property - def mode(self): - return self.meta.get('mode') - - @property - def polarisation(self): - return self.meta.get('polarisation') - - @classmethod - def is_datasource_for(cls, data, meta, **kwargs): - return meta['instrument'] == 'ILOFAR' diff --git a/radiospectra/spectrogram2/spectrogram.py b/radiospectra/spectrogram2/spectrogram.py deleted file mode 100644 index 5dea66a..0000000 --- a/radiospectra/spectrogram2/spectrogram.py +++ /dev/null @@ -1,773 +0,0 @@ -import os -import glob -import gzip -import struct -import pathlib -import warnings -from pathlib import Path -from collections import OrderedDict -from urllib.request import Request, urlopen - -import cdflib -import matplotlib.dates as mdates -import numpy as np -import pandas as pd -from matplotlib import pyplot as plt -from matplotlib.image import NonUniformImage -from scipy.io.idl import readsav - -import astropy.units as u -from astropy.io.fits import Header -from astropy.time import Time -from astropy.visualization import quantity_support -from sunpy import log -from sunpy.data import cache -from sunpy.io import fits -from sunpy.net import attrs as a -from sunpy.time import parse_time -from sunpy.util.datatype_factory_base import ( - BasicRegistrationFactory, - MultipleMatchError, - NoMatchError, - ValidationFunctionError, -) -from sunpy.util.exceptions import SunpyUserWarning -from sunpy.util.functools import seconddispatch -from sunpy.util.metadata import MetaDict -from sunpy.util.util import expand_list - -SUPPORTED_ARRAY_TYPES = (np.ndarray,) -try: - import dask.array - SUPPORTED_ARRAY_TYPES += (dask.array.Array,) -except ImportError: - pass - -quantity_support() - -__all__ = ['SpectrogramFactory', 'Spectrogram', 'GenericSpectrogram'] - - -class NoSpectrogramInFileError(Exception): - pass - - -class SpectraMetaValidationError(AttributeError): - pass - - -# Can use sunpy.until.io once min version has these -def is_url(arg): - try: - urlopen(arg) - except Exception: - return False - return True - - -# In python<3.8 paths with un-representable chars (ie. '*' on windows) -# raise an error, so make our own version that returns False instead of -# erroring. These can be removed when we support python >= 3.8 -# https://docs.python.org/3/library/pathlib.html#methods -def is_file(path): - try: - return path.is_file() - except Exception: - return False - - -def is_dir(path): - try: - return path.is_dir() - except Exception: - return False - - -def possibly_a_path(obj): - """ - Check if ``obj`` can be coerced into a pathlib.Path object. - Does *not* check if the path exists. - """ - try: - pathlib.Path(obj) - return True - except Exception: - return False - - -def parse_path(path, f, **kwargs): - """ - Read in a series of files at *path* using the function *f*. - - Parameters - ---------- - path : pathlib.Path - f : callable - Must return a list of read-in data. - kwargs : - Additional keyword arguments are handed to ``f``. - - Returns - ------- - list - List of files read in by ``f``. - """ - if not isinstance(path, os.PathLike): - raise ValueError('path must be a pathlib.Path object') - path = path.expanduser() - if is_file(path): - read_file = f(path, **kwargs) - if isinstance(read_file, tuple): - read_file = [read_file] - return read_file - elif is_dir(path): - read_files = [] - for afile in sorted(path.glob('*')): - read_files += f(afile, **kwargs) - return read_files - elif glob.glob(str(path)): - read_files = [] - for afile in sorted(glob.glob(str(path))): - read_files += f(afile, **kwargs) - return read_files - else: - raise ValueError(f'Did not find any files at {path}') - - -class PcolormeshPlotMixin: - """ - Class provides plotting functions using `~pcolormesh`. - """ - def plot(self, axes=None, **kwargs): - """ - Plot the spectrogram. - - Parameters - ---------- - axes : `matplotlib.axis.Axes`, optional - The axes where the plot will be added. - kwargs : - Arguments pass to the plot call `pcolormesh`. - - Returns - ------- - """ - if axes is None: - fig, axes = plt.subplots() - else: - fig = axes.get_figure() - - if hasattr(self.data, 'value'): - data = self.data.value - else: - data = self.data - - title = f'{self.observatory}, {self.instrument}' - if self.instrument != self.detector: - title = f'{title}, {self.detector}' - - axes.set_title(title) - axes.plot(self.times.datetime[[0, -1]], self.frequencies[[0, -1]], - linestyle='None', marker='None') - if (self.times.shape[0] == self.data.shape[0] - and self.frequencies.shape[0] == self.data.shape[1]): - axes.pcolormesh(self.times.datetime, self.frequencies.value, data, shading='auto', - **kwargs) - else: - axes.pcolormesh(self.times.datetime, self.frequencies.value, data[:-1, :-1], - shading='auto', **kwargs) - axes.set_xlim(self.times.datetime[0], self.times.datetime[-1]) - locator = mdates.AutoDateLocator(minticks=4, maxticks=8) - formatter = mdates.ConciseDateFormatter(locator) - axes.xaxis.set_major_locator(locator) - axes.xaxis.set_major_formatter(formatter) - fig.autofmt_xdate() - - -class NonUniformImagePlotMixin: - """ - Class provides plotting functions using `NonUniformImage` - """ - def plotim(self, fig=None, axes=None, **kwargs): - if axes is None: - fig, axes = plt.subplots() - - im = NonUniformImage(axes, interpolation='none', **kwargs) - im.set_data(mdates.date2num(self.times.datetime), self.frequencies.value, self.data) - axes.images.append(im) - - -class GenericSpectrogram(PcolormeshPlotMixin, NonUniformImagePlotMixin): - """ - Base spectrogram class all new spectrograms will inherit. - - Attributes - ---------- - meta : `dict-like` - Meta data for the spectrogram. - data : `numpy.ndarray` - The spectrogram data itself a 2D array. - """ - _registry = {} - - def __init_subclass__(cls, **kwargs): - super().__init_subclass__(**kwargs) - if hasattr(cls, 'is_datasource_for'): - cls._registry[cls] = cls.is_datasource_for - - def __init__(self, data, meta, **kwargs): - self.data = data - self.meta = meta - - self._validate_meta() - - @property - def observatory(self): - """ - The name of the observatory which recorded the spectrogram. - """ - return self.meta['observatory'].upper() - - @property - def instrument(self): - """ - The name of the instrument which recorded the spectrogram. - """ - return self.meta['instrument'].upper() - - @property - def detector(self): - """ - The detector which recorded the spectrogram. - """ - return self.meta['detector'].upper() - - @property - def start_time(self): - """ - The start time of the spectrogram. - """ - return self.meta['start_time'] - - @property - def end_time(self): - """ - The end time of the spectrogram. - """ - return self.meta['end_time'] - - @property - def wavelength(self): - """ - The wavelength range of the spectrogram. - """ - return self.meta['wavelength'] - - @property - def times(self): - """ - The times of the spectrogram. - """ - return self.meta['times'] - - @property - def frequencies(self): - """ - The frequencies of the spectrogram. - """ - return self.meta['freqs'] - - def _validate_meta(self): - """ - Validates the meta-information associated with a Spectrogram. - - This method includes very basic validation checks which apply to - all of the kinds of files that radiospectra can read. Datasource-specific - validation should be handled in the relevant file in the - radiospectra.spectrogram2.sources file. - """ - msg = 'Spectrogram coordinate units for {} axis not present in metadata.' - err_message = [] - for i, ax in enumerate(['times', 'freqs']): - if self.meta.get(ax) is None: - err_message.append(msg.format(ax)) - - if err_message: - raise SpectraMetaValidationError('\n'.join(err_message)) - - def __repr__(self): - return (f'<{self.__class__.__name__} {self.observatory}, {self.instrument}, {self.detector}' - f' {self.wavelength.min} - {self.wavelength.max},' - f' {self.start_time.isot} to {self.end_time.isot}>') - - -class SpectrogramFactory(BasicRegistrationFactory): - """ - A factory for generating spectrograms. - - Parameters - ---------- - \\*inputs - `str` or `pathlib.Path` to the file. - - Returns - ------- - `radiospectra.spectrogram2.Spectrogram` - The spectrogram for the give file - """ - def _validate_meta(self, meta): - """ - Validate a meta argument. - """ - if isinstance(meta, Header): - return True - elif isinstance(meta, dict): - return True - else: - return False - - def _parse_args(self, *args, silence_errors=False, **kwargs): - """ - Parses an args list into data-header pairs. - args can contain any mixture of the following entries: - * tuples of data,header - * data, header not in a tuple - * data, wcs object in a tuple - * data, wcs object not in a tuple - * filename, as a str or pathlib.Path, which will be read - * directory, as a str or pathlib.Path, from which all files will be read - * glob, from which all files will be read - * url, which will be downloaded and read - * lists containing any of the above. - - Example - ------- - self._parse_args(data, header, - (data, header), - ['file1', 'file2', 'file3'], - 'file4', - 'directory1', - '*.fits') - """ - # Account for nested lists of items - args = expand_list(args) - - # Sanitise the input so that each 'type' of input corresponds to a different - # class, so single dispatch can be used later - nargs = len(args) - i = 0 - while i < nargs: - arg = args[i] - if isinstance(arg, SUPPORTED_ARRAY_TYPES): - # The next two items are data and a header - data = args.pop(i) - header = args.pop(i) - args.insert(i, (data, header)) - nargs -= 1 - elif isinstance(arg, str) and is_url(arg): - # Replace URL string with a Request object to dispatch on later - args[i] = Request(arg) - elif possibly_a_path(arg): - # Replace path strings with Path objects - args[i] = pathlib.Path(arg) - i += 1 - - # Parse the arguments - # Note that this list can also contain GenericSpectrogram if they are directly given to - # the factory - data_header_pairs = [] - for arg in args: - try: - data_header_pairs += self._parse_arg(arg, **kwargs) - except NoSpectrogramInFileError as e: - if not silence_errors: - raise - warnings.warn(f"One of the arguments failed to parse with error: {e}", - SunpyUserWarning) - - return data_header_pairs - - # Note that post python 3.8 this can be @functools.singledispatchmethod - @seconddispatch - def _parse_arg(self, arg, **kwargs): - """ - Take a factory input and parse into (data, header) pairs. - Must return a list, even if only one pair is returned. - """ - raise ValueError(f"Invalid input: {arg}") - - @_parse_arg.register(tuple) - def _parse_tuple(self, arg, **kwargs): - # Data-header - data, header = arg - pair = data, header - if self._validate_meta(header): - pair = (data, OrderedDict(header)) - return [pair] - - @_parse_arg.register(GenericSpectrogram) - def _parse_map(self, arg, **kwargs): - return [arg] - - @_parse_arg.register(Request) - def _parse_url(self, arg, **kwargs): - url = arg.full_url - path = str(cache.download(url).absolute()) - pairs = self._read_file(path, **kwargs) - return pairs - - @_parse_arg.register(pathlib.Path) - def _parse_path(self, arg, **kwargs): - return parse_path(arg, self._read_file, **kwargs) - - def __call__(self, *args, silence_errors=False, **kwargs): - """ - Method for running the factory. Takes arbitrary arguments and - keyword arguments and passes them to a sequence of pre-registered types - to determine which is the correct spectrogram-type to build. - Arguments args and kwargs are passed through to the validation - function and to the constructor for the final type. For spectrogram types, - validation function must take a data-header pair as an argument. - - Parameters - ---------- - silence_errors : `bool`, optional - If set, ignore data-header pairs which cause an exception. - Default is ``False``. - - Notes - ----- - Extra keyword arguments are passed through to `sunpy.io.read_file` such - as `memmap` for FITS files. - """ - data_header_pairs = self._parse_args(*args, silence_errors=silence_errors, **kwargs) - new_maps = list() - - # Loop over each registered type and check to see if WidgetType - # matches the arguments. If it does, use that type. - for pair in data_header_pairs: - if isinstance(pair, GenericSpectrogram): - new_maps.append(pair) - continue - data, header = pair - meta = MetaDict(header) - - try: - new_map = self._check_registered_widgets(data, meta, **kwargs) - new_maps.append(new_map) - except (NoMatchError, MultipleMatchError, - ValidationFunctionError, SpectraMetaValidationError) as e: - if not silence_errors: - raise - warnings.warn(f'One of the data, header pairs failed to validate with: {e}', - SunpyUserWarning) - - if not len(new_maps): - raise RuntimeError('No maps loaded') - - if len(new_maps) == 1: - return new_maps[0] - - return new_maps - - def _check_registered_widgets(self, data, meta, **kwargs): - - candidate_widget_types = list() - - for key in self.registry: - - # Call the registered validation function for each registered class - if self.registry[key](data, meta, **kwargs): - candidate_widget_types.append(key) - - n_matches = len(candidate_widget_types) - - if n_matches == 0: - if self.default_widget_type is None: - raise NoMatchError("No types match specified arguments and no default is set.") - else: - candidate_widget_types = [self.default_widget_type] - elif n_matches > 1: - raise MultipleMatchError("Too many candidate types identified " - f"({candidate_widget_types}). " - "Specify enough keywords to guarantee unique type " - "identification.") - - # Only one is found - WidgetType = candidate_widget_types[0] - - return WidgetType(data, meta, **kwargs) - - def _read_file(self, file, **kwargs): - file = Path(file) - extensions = file.suffixes - first_extension = extensions[0].lower() - if first_extension == '.dat': - return self._read_dat(file) - elif first_extension in ('.r1', '.r2'): - return self._read_idl_sav(file, instrument='waves') - elif first_extension == '.cdf': - return self._read_cdf(file) - elif first_extension == '.srs': - return self._read_srs(file) - elif first_extension in ('.fits', '.fit', '.fts', 'fit.gz'): - return self._read_fits(file) - else: - raise ValueError(f'Extension {extensions[0]} not supported.') - - @staticmethod - def _read_dat(file): - if 'swaves' in file.name: - name, prod, date, spacecraft, receiver = file.stem.split('_') - # frequency range - freqs = np.genfromtxt(file, max_rows=1) * u.kHz - # bg which is already subtracted from data - bg = np.genfromtxt(file, skip_header=1, max_rows=1) - # data - data = np.genfromtxt(file, skip_header=2) - times = data[:, 0] * u.min - data = data[:, 1:].T - - meta = {'instrument': name, 'observatory': f'STEREO {spacecraft.upper()}', - 'product': prod, 'start_time': Time.strptime(date, '%Y%m%d'), - 'wavelength': a.Wavelength(freqs[0], freqs[-1]), 'detector': receiver, - 'freqs': freqs, 'background': bg} - - meta['times'] = meta['start_time'] + times - meta['end_time'] = meta['start_time'] + times[-1] - - return data, meta - elif 'bst' in file.name: - def subband_to_freq(subband, obs_mode): - """ - Converts LOFAR single station subbands to frequency - - Parameters - ---------- - subband : `int` - Subband number - obs_mode : `int` - Observation mode 3, 5, 7 - - Return - ------ - freq - frequency as astropy.units.Quantity (MHz) - """ - nyq_zone_dict = {3: 1, 5: 2, 7: 3} - nyq_zone = nyq_zone_dict[obs_mode] - clock_dict = {3: 200, 4: 160, 5: 200, 6: 160, 7: 200} # MHz - clock = clock_dict[obs_mode] - freq = (nyq_zone - 1 + subband / 512) * (clock / 2) - return freq * u.MHz # MHz - - subbands = (np.arange(54, 454, 2), np.arange(54, 454, 2), np.arange(54, 230, 2)) - num_subbands = 488 - - data = np.fromfile(file) - polarisation = file.stem[-1] - - num_times = data.shape[0] / num_subbands - if not num_times.is_integer(): - log.warning('BST file seems incomplete dropping incomplete frequencies') - num_times = np.floor(num_times).astype(int) - truncate = num_times * num_subbands - data = data[:truncate] - data = data.reshape(-1, num_subbands).T # (Freq x Time).T = (Time x Freq) - dt = np.arange(num_times) * 1 * u.s - start_time = Time.strptime(file.name.split('_bst')[0], "%Y%m%d_%H%M%S") - times = start_time + dt - - obs_mode = (3, 5, 7) - - freqs = [subband_to_freq(sb, mode) for sb, mode in zip(subbands, obs_mode)] - - # 1st 200 sbs mode 3, next 200 sbs mode 5, last 88 sbs mode 7 - spec = {0: data[:200, :], 1: data[200:400, :], 2: data[400:, :]} - data_header_pairs = [] - for i in range(3): - meta = {'instrument': 'ILOFAR', 'observatory': 'Birr (IE613)', - 'start_time': times[0], 'mode': obs_mode[i], - 'wavelength': a.Wavelength(freqs[i][0], freqs[i][-1]), - 'freqs': freqs[i], 'times': times, 'end_time': times[-1], - 'detector': 'ILOFAR', 'polarisation': polarisation} - - data_header_pairs.append((spec[i], meta)) - - return data_header_pairs - - @staticmethod - def _read_srs(file): - - with file.open('rb') as buff: - data = buff.read() - if file.suffixes[-1] == '.gz': - data = gzip.decompress(data) - - # Data is store as a series of records made of different numbers of bytes - # General header information - # 1 Year (last 2 digits) Byte integer (unsigned) - # 2 Month number (1 to 12) " - # 3 Day (1 to 31) " - # 4 Hour (0 to 23 UT) " - # 5 Minute (0 to 59) " - # 6 Second at start of scan (0 to 59) " - # 7 Site Number (0 to 255) " - # 8 Number of bands in the record (2) " - # - # Band 1 (A-band) header information - # 9,10 Start Frequency (MHz) Word integer (16 bits) - # 11,12 End Frequency (MHz) " - # 13,14 Number of bytes in data record (401)" - # 15 Analyser reference level Byte integer - # 16 Analyser attenuation (dB) " - # - # Band 2 (B-band) header information - # 17-24 As for band 1 - # - # Spectrum Analyser data - # 25-425 401 data bytes for band 1 (A-band) - # 426-826 401 data bytes for band 2 (B-band) - record_struc = struct.Struct('B'*8 + 'H'*3 + 'B'*2 - + 'H'*3 + 'B'*2 + 'B' * 401 + 'B' * 401) - - records = record_struc.iter_unpack(data) - - # Map of numeric records to locations - site_map = {1: 'Palehua', 2: 'Holloman', 3: 'Learmonth', 4: 'San Vito'} - - df = pd.DataFrame([(*r[:18], np.array(r[18:419]), np.array(r[419:820])) - for r in records]) - df.columns = ['year', 'month', 'day', 'hour', 'minute', 'second', 'site', ' num_bands', - 'start_freq1', 'end_freq1', 'num_bytes1', 'analyser_ref1', - 'analyser_atten1', 'start_freq2', 'end_freq2', 'num_bytes2', - 'analyser_ref2', 'analyser_atten2', 'spec1', 'spec2'] - - # Hack to make to_datetime work - earliest dates seem to be 2000 and won't be - # around in 3000! - df['year'] = df['year'] + 2000 - df['time'] = pd.to_datetime(df[['year', 'month', 'day', 'hour', 'minute', 'second']]) - - # Equations taken from document - n = np.arange(1, 402) - freq_a = (25 + 50 * (n - 1) / 400) * u.MHz - freq_b = (75 + 105 * (n - 1) / 400) * u.MHz - freqs = np.hstack([freq_a, freq_b]) - - data = np.hstack([np.vstack(df[name].to_numpy()) for name in ['spec1', 'spec2']]).T - times = Time(df['time']) - - meta = {'instrument': 'RSTN', 'observatory': site_map[df['site'][0]], - 'start_time': times[0], 'end_time': times[-1], 'detector': 'RSTN', - 'wavelength': a.Wavelength(freqs[0], freqs[-1]), 'freqs': freqs, 'times': times} - - return data, meta - - @staticmethod - def _read_cdf(file): - cdf = cdflib.CDF(file) - cdf_meta = cdf.globalattsget() - if (cdf_meta.get('Project', '') == 'PSP' - and cdf_meta.get('Source_name') == 'PSP_FLD>Parker Solar Probe FIELDS' - and 'Radio Frequency Spectrometer' in cdf_meta.get('Descriptor')): - short, _long = cdf_meta['Descriptor'].split('>') - detector = short[4:].lower() - - times, data, freqs = [cdf.varget(name) for name in - [f'epoch_{detector}_auto_averages_ch0_V1V2', - f'psp_fld_l2_rfs_{detector}_auto_averages_ch0_V1V2', - f'frequency_{detector}_auto_averages_ch0_V1V2']] - - times = Time(times << u.ns, format='cdf_tt2000') - freqs = freqs[0, :] << u.Hz - data = data.T << u.Unit('Volt**2/Hz') - - meta = { - 'cdf_meta': cdf_meta, - 'detector': detector, - 'instrument': 'FIELDS/RFS', - 'observatory': 'PSP', - 'start_time': times[0], - 'end_time': times[-1], - 'wavelength': a.Wavelength(freqs.min(), freqs.max()), - 'times': times, - 'freqs': freqs - } - return data, meta - - @staticmethod - def _read_fits(file): - hd_pairs = fits.read(file) - - if 'e-CALLISTO' in hd_pairs[0].header.get('CONTENT', ''): - data = hd_pairs[0].data - times = hd_pairs[1].data['TIME'].flatten() * u.s - freqs = hd_pairs[1].data['FREQUENCY'].flatten() * u.MHz - start_time = parse_time(hd_pairs[0].header['DATE-OBS'] - + ' ' + hd_pairs[0].header['TIME-OBS']) - end_time = parse_time(hd_pairs[0].header['DATE-END'] - + ' ' + hd_pairs[0].header['TIME-END']) - times = start_time + times - meta = { - 'fits_meta': hd_pairs[0].header, - 'detector': 'e-CALLISTO', - 'instrument': 'e-CALLISTO', - 'observatory': hd_pairs[0].header['INSTRUME'], - 'start_time': start_time, - 'end_time': end_time, - 'wavelength': a.Wavelength(freqs.min(), freqs.max()), - 'times': times, - 'freqs': freqs - } - return data, meta - elif hd_pairs[0].header.get('TELESCOP', '') == 'EOVSA': - times = Time(hd_pairs[2].data['mjd'] + hd_pairs[2].data['time'] / 1000.0 / 86400., - format='mjd') - freqs = hd_pairs[1].data['sfreq'] * u.GHz - data = hd_pairs[0].data - start_time = parse_time(hd_pairs[0].header['DATE_OBS']) - end_time = parse_time(hd_pairs[0].header['DATE_END']) - - meta = { - 'fits_meta': hd_pairs[0].header, - 'detector': 'EOVSA', - 'instrument': 'EOVSA', - 'observatory': 'Owens Valley', - 'start_time': start_time, - 'end_time': end_time, - 'wavelength': a.Wavelength(freqs.min(), freqs.max()), - 'times': times, - 'freqs': freqs - } - return data, meta - - @staticmethod - def _read_idl_sav(file, instrument=None): - data = readsav(file) - if instrument == 'waves': - # See https://solar-radio.gsfc.nasa.gov/wind/one_minute_doc.html - data_array = data['arrayb'] - # frequency range - if file.suffix == '.R1': - freqs = np.linspace(20, 1040, 256) * u.kHz - receiver = 'RAD1' - elif file.suffix == '.R2': - freqs = np.linspace(1.075, 13.825, 256) * u.MHz - receiver = 'RAD2' - # bg which is already subtracted from data ? - bg = data_array[:, -1] - data = data_array[:, :-1] - - start_time = Time.strptime(file.stem, '%Y%m%d') - end_time = start_time + 86399*u.s - times = start_time + (np.arange(1440) * 60 + 30)*u.s - - meta = {'instrument': 'WAVES', 'observatory': 'WIND', 'start_time': start_time, - 'end_time': end_time, 'wavelength': a.Wavelength(freqs[0], freqs[-1]), - 'detector': receiver, 'freqs': freqs, 'times': times, 'background': bg} - - return data, meta - - -Spectrogram = SpectrogramFactory(registry=GenericSpectrogram._registry, - default_widget_type=GenericSpectrogram) diff --git a/radiospectra/spectrogram2/tests/test_ilofar.py b/radiospectra/spectrogram2/tests/test_ilofar.py deleted file mode 100644 index dd1fa0e..0000000 --- a/radiospectra/spectrogram2/tests/test_ilofar.py +++ /dev/null @@ -1,28 +0,0 @@ -from pathlib import Path -from unittest import mock - -import numpy as np - -from radiospectra.spectrogram2.spectrogram import Spectrogram - - -@mock.patch('numpy.fromfile') -@mock.patch('radiospectra.spectrogram2.spectrogram.is_file') -def test_ilofar(mock_is_file, mock_fromfile): - mock_fromfile.return_value = np.ones(10117216) - mock_is_file.return_value = True - - spec = Spectrogram(Path('20180602_063247_bst_00X.dat')) - assert len(spec) == 3 - assert spec[0].polarisation == 'X' - assert spec[0].start_time.iso == '2018-06-02 06:32:47.000' - assert spec[0].end_time.iso == '2018-06-02 12:18:18.000' - assert spec[0].mode == 3 - assert spec[0].frequencies[0].to_value('MHz') == 10.546875 - assert spec[0].frequencies[-1].to_value('MHz') == 88.28125 - assert spec[1].mode == 5 - assert spec[1].frequencies[0].to_value('MHz') == 110.546875 - assert spec[1].frequencies[-1].to_value('MHz') == 188.28125 - assert spec[2].mode == 7 - assert spec[2].frequencies[0].to_value('MHz') == 210.546875 - assert spec[2].frequencies[-1].to_value('MHz') == 244.53125 diff --git a/radiospectra/utils.py b/radiospectra/utils.py new file mode 100644 index 0000000..6f356ac --- /dev/null +++ b/radiospectra/utils.py @@ -0,0 +1,29 @@ +import astropy.units as u + +__all__ = ["subband_to_freq"] + + +def subband_to_freq(subband, obs_mode): + """ + Converts LOFAR single station subbands to frequency + + Parameters + ---------- + subband : `int` + Subband number. + obs_mode : `int` + Observation mode 3, 5, 7. + + Return + ------ + `astropy.units.Quantity` + Frequency in MHz + """ + nyq_zone_dict = {3: 1, 5: 2, 7: 3} + if obs_mode not in nyq_zone_dict: + raise ValueError(f"Observation mode {obs_mode} not supported, only 3, 5, 7 are supported.") + nyq_zone = nyq_zone_dict[obs_mode] + clock_dict = {3: 200, 4: 160, 5: 200, 6: 160, 7: 200} # MHz + clock = clock_dict[obs_mode] + freq = (nyq_zone - 1 + subband / 512) * (clock / 2) + return freq * u.MHz # MHz From fe59bd6e0d8aa01e09f0098ff8f698d44cca952f Mon Sep 17 00:00:00 2001 From: Nabil Freij Date: Fri, 17 Nov 2023 22:19:31 -0800 Subject: [PATCH 8/8] Fix bad merge and bug(?) --- radiospectra/net/__init__.py | 2 +- radiospectra/net/sources/ilofar.py | 28 ++++++++++--------- radiospectra/spectrogram/sources/__init__.py | 2 ++ .../spectrogram/spectrogram_factory.py | 3 +- 4 files changed, 19 insertions(+), 16 deletions(-) diff --git a/radiospectra/net/__init__.py b/radiospectra/net/__init__.py index 778cf68..a71ea41 100644 --- a/radiospectra/net/__init__.py +++ b/radiospectra/net/__init__.py @@ -8,7 +8,7 @@ from radiospectra.net.sources.wind import WAVESClient __all__ = [ - "CALLISTOClient", + "eCALLISTOClient", "EOVSAClient", "RFSClient", "SWAVESClient", diff --git a/radiospectra/net/sources/ilofar.py b/radiospectra/net/sources/ilofar.py index a0aa5d8..c299fce 100644 --- a/radiospectra/net/sources/ilofar.py +++ b/radiospectra/net/sources/ilofar.py @@ -26,21 +26,23 @@ class ILOFARMode357Client(GenericClient): >>> results = Fido.search(a.Time("2021/09/01", "2021/09/21"), ... a.Instrument('ILOFAR')) # doctest: +REMOTE_DATA >>> results #doctest: +REMOTE_DATA - Results from 1 Provider: - Start Time End Time ... PolType Polarisation - ----------------------- ----------------------- ... ------- ------------ - 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ... X X - 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ... X Y - 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ... X X - 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ... X Y - 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ... X X - 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ... X Y - 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ... X X - 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ... X Y - 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ... X X - 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ... X Y + 10 Results from the ILOFARMode357Client: + + Start Time End Time ... Provider Polarisation + ----------------------- ----------------------- ... -------- ------------ + 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ... ILOFAR X + 2021-09-14 07:39:13.000 2021-09-14 07:39:13.999 ... ILOFAR Y + 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ... ILOFAR X + 2021-09-01 08:07:29.000 2021-09-01 08:07:29.999 ... ILOFAR Y + 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ... ILOFAR X + 2021-09-07 08:07:52.000 2021-09-07 08:07:52.999 ... ILOFAR Y + 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ... ILOFAR X + 2021-09-08 08:04:07.000 2021-09-08 08:04:07.999 ... ILOFAR Y + 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ... ILOFAR X + 2021-09-08 10:34:31.000 2021-09-08 10:34:31.999 ... ILOFAR Y """ diff --git a/radiospectra/spectrogram/sources/__init__.py b/radiospectra/spectrogram/sources/__init__.py index 09316f4..60336ee 100644 --- a/radiospectra/spectrogram/sources/__init__.py +++ b/radiospectra/spectrogram/sources/__init__.py @@ -8,6 +8,7 @@ from ..spectrogram_factory import Spectrogram # NOQA from .callisto import * # NOQA from .eovsa import * # NOQA +from .ilofar357 import * # NOQA from .psp_rfs import * # NOQA from .rpw import * # NOQA from .rstn import * # NOQA @@ -21,4 +22,5 @@ "EOVSASpectrogram", "RSTNSpectrogram", "RPWSpectrogram", + "ILOFARMode357Spectrogram", ] diff --git a/radiospectra/spectrogram/spectrogram_factory.py b/radiospectra/spectrogram/spectrogram_factory.py index 06b411c..5988288 100644 --- a/radiospectra/spectrogram/spectrogram_factory.py +++ b/radiospectra/spectrogram/spectrogram_factory.py @@ -238,7 +238,7 @@ def _read_file(self, file, **kwargs): extensions = file.suffixes first_extension = extensions[0].lower() if first_extension == ".dat": - return [self._read_dat(file)] + return self._read_dat(file) elif first_extension in (".r1", ".r2"): return [self._read_idl_sav(file, instrument="waves")] elif first_extension == ".cdf": @@ -318,7 +318,6 @@ def _read_dat(file): } data_header_pairs.append((spec[i], meta)) - return data_header_pairs else: raise ValueError(f"File {file} not supported.")