diff --git a/testfile cartopy.py b/testfile cartopy.py new file mode 100644 index 0000000..75627c1 --- /dev/null +++ b/testfile cartopy.py @@ -0,0 +1,186 @@ +#Lets import some stuff! +import boto3 +from botocore import UNSIGNED +from botocore.client import Config +from datetime import timedelta, datetime, timezone +import os +import pyart +from matplotlib import pyplot as plt +import tempfile +import numpy as np +import io +from PIL import Image + +import cartopy.crs as ccrs +from cartopy.io.img_tiles import OSM + +def printdebug(msg): + print(f'DEBUG: {msg}') + +#Helper function for the search +def _nearestDate(dates, pivot): + return min(dates, key=lambda x: abs(x - pivot)) + +def _getRadarFileNames(site, datetime_t): + return datetimes + +def get_radar_frames(site, datetime_t, num_frames): + #First create the query string for the bucket knowing + #how NOAA and AWS store the data + + currentDateString = datetime_t.strftime('%Y/%m/%d/') + site + previousDateString = (datetime_t - timedelta(days=1)).strftime('%Y/%m/%d/') + site + + conn = boto3.resource('s3', config=Config(signature_version=UNSIGNED)) + bucket = conn.Bucket('noaa-nexrad-level2') + + bucketListCurrentDate = list(bucket.objects.filter(Prefix = currentDateString)) + bucketListPreviousDate = list(bucket.objects.filter(Prefix = previousDateString)) + + #we are going to create a list of keys and datetimes to allow easy searching + + keys = [] + datetimes = [] + + #populate the list + printdebug('Getting nearby dates') + for i in range(len(bucketListCurrentDate)): + this_str = str(bucketListCurrentDate[i].key) + if 'gz' in this_str: + endme = this_str[-22:-4] + fmt = '%Y%m%d_%H%M%S_V0' + dt = datetime.strptime(endme, fmt).replace(tzinfo=timezone.utc) + datetimes.append(dt) + keys.append(bucketListCurrentDate[i]) + + if this_str[-3::] == 'V06': + endme = this_str[-19::] + fmt = '%Y%m%d_%H%M%S_V06' + dt = datetime.strptime(endme, fmt).replace(tzinfo=timezone.utc) + datetimes.append(dt) + keys.append(bucketListCurrentDate[i]) + + #datetimes = _getRadarFileNames(site, datetime_t) + printdebug('Getting nearest datetime') + closest_datetime = _nearestDate(datetimes, datetime_t) + index = datetimes.index(closest_datetime) + + radarFrames = [0] * num_frames + i = index + printdebug('Pulling all the radar images requested') + while i > (index - num_frames) and i >= 0: + radarFrames[num_frames - ((index - i) + 1)] = get_radar_from_aws(keys[i].Object()) + i = i - 1 + + return radarFrames + +def get_radar_from_aws(aws_object): + """ + Get the closest volume of NEXRAD data to a particular datetime. + Parameters + ---------- + site : string + four letter radar designation + datetime_t : datetime + desired date time + Returns + ------- + radar : Py-ART Radar Object + Radar closest to the queried datetime + """ + + #find the closest available radar to your datetime + + #closest_datetime = _nearestDate(datetimes, datetime_t) + #index = datetimes.index(closest_datetime) + + #localfile = tempfile.NamedTemporaryFile() + #keys[index].Object().download_file(localfile.name) + #aws_object.download_file(localfile.name) + #radar = pyart.io.read(localfile.name) + filename = aws_object.key.replace('/', '_') + if not os.path.isfile(f'cache/rawradar/{filename}'): + printdebug('Cache Miss') + aws_object.download_file(f'cache/rawradar/{filename}') + radar = pyart.io.read(f'cache/rawradar/{filename}') + return radar + +def convert_fig_to_img(fig): + fig.canvas.draw() + (w,h) = fig.canvas.get_width_height() + img = Image.frombytes('RGBA', fig.canvas.get_width_height(), np.frombuffer(fig.canvas.buffer_rgba(), dtype=np.uint8).reshape((h,w,4))) + return img + +def generate_radar_image(radar_data, tiler, mercator): + #determine max_lat / max_lon etc from position of radar + max_lat = radar_data.latitude["data"][0] + 4.5 + min_lat = radar_data.latitude["data"][0] - 4.5 + max_lon = radar_data.longitude["data"][0] + 5.5 + min_lon = radar_data.longitude["data"][0] - 5.5 + + # Used to be .5 + lal = np.arange(min_lat, max_lat, .2) + lol = np.arange(min_lon, max_lon, .1) + + printdebug('Creating RadarMapDisplay') + display = pyart.graph.RadarMapDisplay(radar_data) + + projection = ccrs.LambertConformal( + central_latitude=radar_data.latitude["data"][0], + central_longitude=radar_data.longitude["data"][0], + ) + + printdebug('Processing data to clean it up a bit') + + radar_data.fields['reflectivity']['data'][:, -10:] = np.ma.masked + + gatefilter = pyart.filters.GateFilter(radar_data) + gatefilter.exclude_transition() + gatefilter.exclude_masked("reflectivity") + gatefilter.exclude_below('reflectivity', 16) + + printdebug('Despeckling') + despeckle_gatefilter = pyart.correct.despeckle_field(radar_data, 'reflectivity', gatefilter=gatefilter, size=20) + + fig = plt.figure(figsize = [10,8]) + fig.patch.set_facecolor('none') + ax = fig.add_subplot(1, 1, 1, projection=mercator) + #ax = fig.add_subplot(1, 1, 1) + ax.set_extent([min_lon, max_lon, min_lat, max_lat], crs=ccrs.PlateCarree()) + ax.add_image(tiler, 6) + printdebug('Plotting!') + display.plot_ppi_map('reflectivity', sweep = 0, resolution='10m', + vmin = 10, vmax = 64, + cmap = "ChaseSpectral", + min_lat = min_lat, min_lon = min_lon, + max_lat = max_lat, max_lon = max_lon, + lat_lines = lal, lon_lines = lol, + projection=projection, fig=fig, ax=ax, title_flag=False, + lat_0=radar_data.latitude["data"][0], lon_0=radar_data.longitude["data"][0], + gatefilter=despeckle_gatefilter, add_grid_lines=False, alpha=0.8) + + display.plot_range_ring(radar_data.range["data"][-1] / 1000.0) + + return fig + +fmt = '%Y%m%d_%H%M%S' +currentDatetime = datetime.now(timezone.utc) +num_frames = 10 +secondary_frames = [0] * (num_frames - 1) + +radar = get_radar_frames('KVNX', currentDatetime, num_frames=num_frames) + +tiler = OSM(cache=True) +mercator = tiler.crs + +firstFrameFig = generate_radar_image(radar[0], tiler, mercator) +firstFrame = convert_fig_to_img(firstFrameFig) +i = 0 +for i in range(num_frames - 1): + fig = generate_radar_image(radar[i + 1], tiler, mercator) + secondary_frames[i] = convert_fig_to_img(fig) + # filename = aws_object.key.replace('/', '_') + +datetimeString = datetime.strftime(currentDatetime, fmt) + +firstFrame.save(f'testradar_{datetimeString}.gif', append_images=secondary_frames, duration=200)