Merge 533e2a2 into c19ce9f

Author: pabloitu

csep/__init__.py

@@ -61,7 +61,9 @@
     '__version__'
 ]
 
-def load_stochastic_event_sets(filename, type='csv', format='native', **kwargs):
+
+def load_stochastic_event_sets(filename, type='csv', format='native',
+                               **kwargs):
     """ General function to load stochastic event sets
 
     This function returns a generator to iterate through a collection of catalogs.
@@ -106,7 +108,8 @@ def load_stochastic_event_sets(filename, type='csv', format='native', **kwargs):
             raise ValueError('format must be either "native" or "csep!')
 
 
-def load_catalog(filename, type='csep-csv', format='native', loader=None, apply_filters=False, **kwargs):
+def load_catalog(filename, type='csep-csv', format='native', loader=None,
+                 apply_filters=False, **kwargs):
     """ General function to load single catalog
 
     See corresponding class documentation for additional parameters.
@@ -122,9 +125,12 @@ def load_catalog(filename, type='csep-csv', format='native', loader=None, apply_
     Returns (:class:`~csep.core.catalogs.AbstractBaseCatalog`)
     """
 
-
-    if type not in ('ucerf3', 'csep-csv', 'zmap', 'jma-csv', 'ingv_horus', 'ingv_emrcmt', 'ndk') and loader is None:
-        raise ValueError("type must be one of the following: ('ucerf3', 'csep-csv', 'zmap', 'jma-csv', 'ndk', 'ingv_horus', 'ingv_emrcmt').")
+    if type not in (
+            'ucerf3', 'csep-csv', 'zmap', 'jma-csv', 'ingv_horus',
+            'ingv_emrcmt',
+            'ndk') and loader is None:
+        raise ValueError(
+            "type must be one of the following: ('ucerf3', 'csep-csv', 'zmap', 'jma-csv', 'ndk', 'ingv_horus', 'ingv_emrcmt').")
 
     # map to correct catalog class, at some point these could be abstracted into configuration file
     # this maps a human readable string to the correct catalog class and the correct loader function
@@ -167,7 +173,8 @@ def load_catalog(filename, type='csep-csv', format='native', loader=None, apply_
         if loader is None:
             loader = class_loader_mapping[type]['loader']
 
-        catalog = catalog_class.load_catalog(filename=filename, loader=loader, **kwargs)
+        catalog = catalog_class.load_catalog(filename=filename, loader=loader,
+                                             **kwargs)
 
     # convert to csep format if needed
     if format == 'native':
@@ -213,12 +220,15 @@ def query_comcat(start_time, end_time, min_magnitude=2.50,
     # Timezone should be in UTC
     t0 = time.time()
     eventlist = readers._query_comcat(start_time=start_time, end_time=end_time,
-                           min_magnitude=min_magnitude,
-                           min_latitude=min_latitude, max_latitude=max_latitude,
-                           min_longitude=min_longitude, max_longitude=max_longitude,
-                           max_depth=max_depth)
+                                      min_magnitude=min_magnitude,
+                                      min_latitude=min_latitude,
+                                      max_latitude=max_latitude,
+                                      min_longitude=min_longitude,
+                                      max_longitude=max_longitude,
+                                      max_depth=max_depth)
     t1 = time.time()
-    comcat = catalogs.CSEPCatalog(data=eventlist, date_accessed=utc_now_datetime(), **kwargs)
+    comcat = catalogs.CSEPCatalog(data=eventlist,
+                                  date_accessed=utc_now_datetime(), **kwargs)
     print("Fetched ComCat catalog in {} seconds.\n".format(t1 - t0))
 
     if apply_filters:
@@ -229,9 +239,13 @@ def query_comcat(start_time, end_time, min_magnitude=2.50,
 
     if verbose:
         print("Downloaded catalog from ComCat with following parameters")
-        print("Start Date: {}\nEnd Date: {}".format(str(comcat.start_time), str(comcat.end_time)))
-        print("Min Latitude: {} and Max Latitude: {}".format(comcat.min_latitude, comcat.max_latitude))
-        print("Min Longitude: {} and Max Longitude: {}".format(comcat.min_longitude, comcat.max_longitude))
+        print("Start Date: {}\nEnd Date: {}".format(str(comcat.start_time),
+                                                    str(comcat.end_time)))
+        print(
+            "Min Latitude: {} and Max Latitude: {}".format(comcat.min_latitude,
+                                                           comcat.max_latitude))
+        print("Min Longitude: {} and Max Longitude: {}".format(
+            comcat.min_longitude, comcat.max_longitude))
         print("Min Magnitude: {}".format(comcat.min_magnitude))
         print(f"Found {comcat.event_count} events in the ComCat catalog.")
 
@@ -266,11 +280,14 @@ def query_bsi(start_time, end_time, min_magnitude=2.50,
     t0 = time.time()
     eventlist = readers._query_bsi(start_time=start_time, end_time=end_time,
                                    min_magnitude=min_magnitude,
-                                   min_latitude=min_latitude, max_latitude=max_latitude,
-                                   min_longitude=min_longitude, max_longitude=max_longitude,
+                                   min_latitude=min_latitude,
+                                   max_latitude=max_latitude,
+                                   min_longitude=min_longitude,
+                                   max_longitude=max_longitude,
                                    max_depth=max_depth)
     t1 = time.time()
-    bsi = catalogs.CSEPCatalog(data=eventlist, date_accessed=utc_now_datetime(), **kwargs)
+    bsi = catalogs.CSEPCatalog(data=eventlist,
+                               date_accessed=utc_now_datetime(), **kwargs)
     print("Fetched BSI catalog in {} seconds.\n".format(t1 - t0))
 
     if apply_filters:
@@ -280,19 +297,25 @@ def query_bsi(start_time, end_time, min_magnitude=2.50,
             bsi = bsi.filter()
 
     if verbose:
-        print("Downloaded catalog from Bollettino Sismico Italiano (BSI) with following parameters")
-        print("Start Date: {}\nEnd Date: {}".format(str(bsi.start_time), str(bsi.end_time)))
-        print("Min Latitude: {} and Max Latitude: {}".format(bsi.min_latitude, bsi.max_latitude))
-        print("Min Longitude: {} and Max Longitude: {}".format(bsi.min_longitude, bsi.max_longitude))
+        print(
+            "Downloaded catalog from Bollettino Sismico Italiano (BSI) with following parameters")
+        print("Start Date: {}\nEnd Date: {}".format(str(bsi.start_time),
+                                                    str(bsi.end_time)))
+        print("Min Latitude: {} and Max Latitude: {}".format(bsi.min_latitude,
+                                                             bsi.max_latitude))
+        print(
+            "Min Longitude: {} and Max Longitude: {}".format(bsi.min_longitude,
+                                                             bsi.max_longitude))
         print("Min Magnitude: {}".format(bsi.min_magnitude))
         print(f"Found {bsi.event_count} events in the BSI catalog.")
 
     return bsi
 
+
 def query_gns(start_time, end_time,  min_magnitude=2.950,
-               min_latitude=-47, max_latitude=-34,
-               min_longitude=164, max_longitude=180,
-               max_depth=45.5,
+              min_latitude=-47, max_latitude=-34,
+              min_longitude=164, max_longitude=180,
+              max_depth=45.5,
               verbose=True,
               apply_filters=False, **kwargs):
     """
@@ -338,6 +361,29 @@ def query_gns(start_time, end_time,  min_magnitude=2.950,
 
     return gns
 
+
+def query_gcmt(start_time, end_time, min_magnitude=5.0,
+               max_depth=None,
+               catalog_id=None,
+               min_latitude=None, max_latitude=None,
+               min_longitude=None, max_longitude=None):
+
+    eventlist = readers._query_gcmt(start_time=start_time,
+                                     end_time=end_time,
+                                     min_magnitude=min_magnitude,
+                                     min_latitude=min_latitude,
+                                     max_latitude=max_latitude,
+                                     min_longitude=min_longitude,
+                                     max_longitude=max_longitude,
+                                     max_depth=max_depth)
+
+    catalog = catalogs.CSEPCatalog(data=eventlist,
+                                   name='gCMT',
+                                   catalog_id=catalog_id,
+                                   date_accessed=utc_now_datetime())
+    return catalog
+
+
 def load_evaluation_result(fname):
     """ Load evaluation result stored as json file
 
@@ -361,7 +407,8 @@ def load_evaluation_result(fname):
             evaluation_type = json_dict['type']
         except:
             evaluation_type = 'default'
-    eval_result = evaluation_result_factory[evaluation_type].from_dict(json_dict)
+    eval_result = evaluation_result_factory[evaluation_type].from_dict(
+        json_dict)
     return eval_result
 
 
@@ -404,15 +451,18 @@ class with the region and magnitude members correctly assigned.
 
     # sanity checks
     if not os.path.exists(fname):
-        raise FileNotFoundError(f"Could not locate file {fname}. Unable to load forecast.")
+        raise FileNotFoundError(
+            f"Could not locate file {fname}. Unable to load forecast.")
     # sanity checks
     if loader is not None and not callable(loader):
-        raise AttributeError("Loader must be callable. Unable to load forecast.")
+        raise AttributeError(
+            "Loader must be callable. Unable to load forecast.")
     extension = os.path.splitext(fname)[-1][1:]
     if extension not in forecast_loader_mapping.keys() and loader is None:
-        raise AttributeError("File extension should be in ('dat','xml','h5','bin') if loader not provided.")
+        raise AttributeError(
+            "File extension should be in ('dat','xml','h5','bin') if loader not provided.")
 
-    if extension in ('xml','h5','bin'):
+    if extension in ('xml', 'h5', 'bin'):
         raise NotImplementedError
 
     # assign default loader
@@ -425,7 +475,8 @@ class with the region and magnitude members correctly assigned.
     return forecast
 
 
-def load_catalog_forecast(fname, catalog_loader=None, format='native', type='ascii', **kwargs):
+def load_catalog_forecast(fname, catalog_loader=None, format='native',
+                          type='ascii', **kwargs):
     """ General function to handle loading catalog forecasts.
 
         Currently, just a simple wrapper, but can contain more complex logic in the future.
@@ -444,10 +495,12 @@ def load_catalog_forecast(fname, catalog_loader=None, format='native', type='asc
     """
     # sanity checks
     if not os.path.exists(fname):
-        raise FileNotFoundError(f"Could not locate file {fname}. Unable to load forecast.")
+        raise FileNotFoundError(
+            f"Could not locate file {fname}. Unable to load forecast.")
     # sanity checks
     if catalog_loader is not None and not callable(catalog_loader):
-        raise AttributeError("Loader must be callable. Unable to load forecast.")
+        raise AttributeError(
+            "Loader must be callable. Unable to load forecast.")
     # factory methods for loading different types of catalogs
     catalog_loader_mapping = {
         'ascii': catalogs.CSEPCatalog.load_ascii_catalogs,
@@ -456,17 +509,18 @@ def load_catalog_forecast(fname, catalog_loader=None, format='native', type='asc
     if catalog_loader is None:
         catalog_loader = catalog_loader_mapping[type]
     # try and parse information from filename and send to forecast constructor
-    if format == 'native' and type=='ascii':
+    if format == 'native' and type == 'ascii':
         try:
             basename = str(os.path.basename(fname.rstrip('/')).split('.')[0])
             split_fname = basename.split('_')
             name = split_fname[0]
-            start_time = strptime_to_utc_datetime(split_fname[1], format="%Y-%m-%dT%H-%M-%S-%f")
+            start_time = strptime_to_utc_datetime(split_fname[1],
+                                                  format="%Y-%m-%dT%H-%M-%S-%f")
             # update kwargs
             _ = kwargs.setdefault('name', name)
             _ = kwargs.setdefault('start_time', start_time)
         except:
             pass
     # create observed_catalog forecast
-    return CatalogForecast(filename=fname, loader=catalog_loader, catalog_format=format, catalog_type=type, **kwargs)
-
+    return CatalogForecast(filename=fname, loader=catalog_loader,
+                           catalog_format=format, catalog_type=type, **kwargs)

csep/utils/iris.py

@@ -0,0 +1,132 @@
+# python imports
+from datetime import datetime
+from urllib import request
+from urllib.parse import urlencode
+
+# PyCSEP imports
+from csep.utils.time_utils import datetime_to_utc_epoch
+
+HOST_CATALOG = "https://service.iris.edu/fdsnws/event/1/query?"
+TIMEOUT = 180
+
+
+def gcmt_search(format='text',
+                starttime=None,
+                endtime=None,
+                updatedafter=None,
+                minlatitude=None,
+                maxlatitude=None,
+                minlongitude=None,
+                maxlongitude=None,
+                latitude=None,
+                longitude=None,
+                maxradius=None,
+                catalog='GCMT',
+                contributor=None,
+                maxdepth=1000,
+                maxmagnitude=10.0,
+                mindepth=-100,
+                minmagnitude=0,
+                offset=1,
+                orderby='time-asc',
+                host=None,
+                verbose=False):
+    """Search the IRIS database for events matching input criteria.
+    This search function is a wrapper around the ComCat Web API described here:
+    https://service.iris.edu/fdsnws/event/1/
+
+    This function returns a list of SummaryEvent objects, described elsewhere in this package.
+    Args:
+        starttime (datetime):
+            Python datetime - Limit to events on or after the specified start time.
+        endtime (datetime):
+            Python datetime - Limit to events on or before the specified end time.
+        updatedafter (datetime):
+           Python datetime - Limit to events updated after the specified time.
+        minlatitude (float):
+            Limit to events with a latitude larger than the specified minimum.
+        maxlatitude (float):
+            Limit to events with a latitude smaller than the specified maximum.
+        minlongitude (float):
+            Limit to events with a longitude larger than the specified minimum.
+        maxlongitude (float):
+            Limit to events with a longitude smaller than the specified maximum.
+        latitude (float):
+            Specify the latitude to be used for a radius search.
+        longitude (float):
+            Specify the longitude to be used for a radius search.
+        maxradius (float):
+            Limit to events within the specified maximum number of degrees
+            from the geographic point defined by the latitude and longitude parameters.
+        catalog (str):
+            Limit to events from a specified catalog.
+        contributor (str):
+            Limit to events contributed by a specified contributor.
+        maxdepth (float):
+            Limit to events with depth less than the specified maximum.
+        maxmagnitude (float):
+            Limit to events with a magnitude smaller than the specified maximum.
+        mindepth (float):
+            Limit to events with depth more than the specified minimum.
+        minmagnitude (float):
+            Limit to events with a magnitude larger than the specified minimum.
+        offset (int):
+            Return results starting at the event count specified, starting at 1.
+        orderby (str):
+            Order the results. The allowed values are:
+            - time order by origin descending time
+            - time-asc order by origin ascending time
+            - magnitude order by descending magnitude
+            - magnitude-asc order by ascending magnitude
+        host (str):
+            Replace default ComCat host (earthquake.usgs.gov) with a custom host.
+    Returns:
+        list: List of SummaryEvent() objects.
+    """
+
+    # getting the inputargs must be the first line of the method!
+    inputargs = locals().copy()
+    newargs = {}
+
+    for key, value in inputargs.items():
+        if value is True:
+            newargs[key] = 'true'
+            continue
+        if value is False:
+            newargs[key] = 'false'
+            continue
+        if value is None:
+            continue
+        newargs[key] = value
+
+    del newargs['verbose']
+
+    events = _search_gcmt(**newargs)
+
+    return events
+
+
+def _search_gcmt(**_newargs):
+    """
+    Performs de-query at ISC API and returns event list and access date
+
+    """
+    paramstr = urlencode(_newargs)
+    url = HOST_CATALOG + paramstr
+    fh = request.urlopen(url, timeout=TIMEOUT)
+    data = fh.read().decode('utf8').split('\n')
+    fh.close()
+    eventlist = []
+    for line in data[1:]:
+        line_ = line.split('|')
+        if len(line_) != 1:
+            id_ = line_[0]
+            time_ = datetime.fromisoformat(line_[1])
+            dt = datetime_to_utc_epoch(time_)
+            lat = float(line_[2])
+            lon = float(line_[3])
+            depth = float(line_[4])
+            mag = float(line_[10])
+            eventlist.append((id_, dt, lat, lon, depth, mag))
+
+    return eventlist