meteortools.rmsutils.multiEventGroundMap

  1# Copyright (C) 2018-2023 Mark McIntyre
  2
  3import datetime
  4import pandas as pd
  5import argparse
  6import os
  7import matplotlib.pyplot as plt
  8from meteortools.utils import greatCircleDistance
  9try:
 10    from wmpl.Utils.PlotMap_OSM import OSMMap
 11except:
 12    print('WMPL not available')
 13
 14RAD2DEG=57.2958
 15
 16
 17def multiEventGroundMap(startdt, enddt, statid=None, shwr=None, outdir=None, minmag=999, fireballsonly=False,
 18                        obslat=None, obslon=None, evtdist=None):
 19    """
 20    Plots a ground track diagram of all events between two dates, with filters by station and shower  
 21
 22    Arguments:  
 23        start:      [string] start date in YYYYMMDD format  
 24        end:        [string] end date in YYYYMMDD format  
 25        statid:     [string] station to filter for. Default all stations.  
 26        shwr:       [string] Filter by shower eg PER. Default All showers.  
 27        outdir:     [string] where to save the file to. if this parameter is omitted, the image will be displayed not saved  
 28        shwr:       [float] Min mag to include.  
 29        fbonly:     [bool] True if only to plot fireballs.
 30        obslat:     [float] Observers latitude: if this is set, obslon and evtdist must be set too
 31        obslon:     [float] Observers longitude
 32        evtdist:    [float] Event distance from observer
 33
 34    Output:  
 35        A jpg map of the detections.   
 36
 37    Note:  
 38        This function reads directly from the UKMON public dataset.  
 39        If obslat is set, then obslon and evtdist must be set, and only events within evtdist of the 
 40        observer will be plotted. 
 41
 42    """
 43    yr = startdt[:4]
 44
 45    cols=['_lat1', '_lng1','_lat2','_lng2','_stream','dtstamp', 'stations', '_amag', 'isfb']
 46    matchfile = f'https://archive.ukmeteors.co.uk/browse/parquet/matches-full-{yr}.parquet.snap'
 47
 48    dta = pd.read_parquet(matchfile, columns=cols)
 49    dta = dta[dta._amag < minmag]
 50    if fireballsonly:
 51        dta = dta[dta.isfb==1]
 52    # filter the data down to just the cols we want
 53    sd = datetime.datetime.strptime(startdt, '%Y%m%d')
 54    sd = sd+datetime.timedelta(hours=12)
 55    if startdt == enddt:
 56        ed = sd + datetime.timedelta(days=1)
 57    else:
 58        ed = datetime.datetime.strptime(enddt, '%Y%m%d')
 59        ed = ed+datetime.timedelta(hours=12)
 60    dta = dta[dta.dtstamp >= sd.timestamp()]
 61    dta = dta[dta.dtstamp <= ed.timestamp()]
 62    
 63    if shwr is not None:
 64        dta = dta[dta._stream==shwr.upper()]
 65    if statid is not None:
 66        dta = dta[dta.stations.str.contains(statid.upper())]
 67    if obslat:
 68        obslat = float(obslat)
 69        obslon = float(obslon)
 70        evtdist = float(evtdist)
 71        dta['evtdist'] = [greatCircleDistance(obslat/RAD2DEG, obslon/RAD2DEG, x/RAD2DEG,y/RAD2DEG) for x,y in zip(dta._lat2, dta._lng2)]
 72        dta = dta[dta.evtdist <= evtdist]
 73
 74    if len(dta) > 1:
 75        plt.clf()
 76        fig = plt.gcf()
 77        fig.set_size_inches(11.6, 8.26)
 78        lat_list = [min(min(dta._lat1), min(dta._lat2))/RAD2DEG, max(max(dta._lat1), max(dta._lat2))/RAD2DEG]
 79        lon_list = [min(min(dta._lng1), min(dta._lng2))/RAD2DEG, max(max(dta._lng1), max(dta._lng2))/RAD2DEG]
 80        # Init the map
 81        #print(lat_list, lon_list)
 82        m = OSMMap(lat_list, lon_list, border_size=50, color_scheme='dark')
 83        lat1s=list(dta._lat1)
 84        lat2s=list(dta._lat2)
 85        lng1s=list(dta._lng1)
 86        lng2s=list(dta._lng2)
 87        for l1, l2, g1, g2 in zip(lat1s, lat2s, lng1s, lng2s):
 88            lats = [l1/RAD2DEG, l2/RAD2DEG]
 89            lons = [g1/RAD2DEG, g2/RAD2DEG]
 90            m.plot(lats, lons, c='r')
 91            m.scatter(l2/RAD2DEG, g2/RAD2DEG, c='k', marker='+', s=50, alpha=0.75)
 92        if outdir is not None:
 93            plt.savefig(os.path.join(outdir, f'{startdt}-{enddt}-{shwr}-{statid}.jpg'))
 94        else:
 95            plt.show()
 96        plt.close()
 97    return 
 98
 99
100if __name__ == '__main__':
101    arg_parser = argparse.ArgumentParser(description="""Plot a ground map of many detections.""",
102        formatter_class=argparse.RawTextHelpFormatter)
103
104    arg_parser.add_argument('start_date', type=str, help='start date in yyyymmdd format')
105    arg_parser.add_argument('end_date', type=str, help='end date in yyyymmdd format')
106    arg_parser.add_argument('-s', '--shower', metavar='SHOWER', type=str,
107        help="Map just this single shower given its code (e.g. PER, ORI, ETA).")
108
109    arg_parser.add_argument('-i', '--stationid', metavar='STATID', help='Station id eg UK0006')
110    arg_parser.add_argument('-o', '--outdir', metavar='OUTDIR', help='Location to save jpg into')
111    arg_parser.add_argument('-m', '--minmag', metavar='MINMAG', type=float, help='Minimum magnitude to filter for')
112    arg_parser.add_argument('-f', '--fbonly', action='store_true', help='Plot only fireballs')
113    arg_parser.add_argument('-l', '--obs_lat', metavar='OBSLAT', help='Observer latitude (degrees')
114    arg_parser.add_argument('-g', '--obs_lon', metavar='OBSLON', help='Observer longitude (degrees)')
115    arg_parser.add_argument('-d', '--event_distance', metavar='EVTDIST', help='Distance from observer (km)')
116
117    cml_args = arg_parser.parse_args()
118    minmag = cml_args.minmag
119    if not minmag:
120        minmag = 999
121    obslat = cml_args.obs_lat
122    if not obslat:
123        obslat = None
124        obslon = None
125        evtdist = None
126    else:
127        if not cml_args.obs_lon or not cml_args.event_distance:
128            print('if providing an observer latitude, must also supply longitude and distance')
129            exit(0)
130        else:
131            obslon = cml_args.obs_lon
132            evtdist = cml_args.event_distance
133    multiEventGroundMap(cml_args.start_date, cml_args.end_date, 
134        cml_args.stationid, cml_args.shower, cml_args.outdir, minmag, cml_args.fbonly, 
135        obslat, obslon, evtdist)