Python MapPlot.drawcities方法代码示例

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def main():
    """Go!"""
    title = 'NOAA MRMS Q3: RADAR + Guage Corrected Rainfall Estimates + NWS Storm Reports'
    mp = MapPlot(sector='custom',
                 north=42.3, east=-93.0, south=41.65, west=-94.1,
                 axisbg='white',
                 titlefontsize=14,
                 title=title,
                 subtitle='Valid: 14 June 2018')

    shp = shapefile.Reader('cities.shp')
    for record in shp.shapeRecords():
        geo = shape(record.shape)
        mp.ax.add_geometries([geo], ccrs.PlateCarree(), zorder=Z_OVERLAY2,
                             facecolor='None', edgecolor='k', lw=2)

    grbs = pygrib.open('MRMS_GaugeCorr_QPE_24H_00.00_20180614-200000.grib2')
    grb = grbs.message(1)
    pcpn = distance(grb['values'], 'MM').value('IN')
    lats, lons = grb.latlons()
    lons -= 360.
    clevs = [0.01, 0.1, 0.3, 0.5, 0.75, 1, 1.5, 2, 2.5, 3, 3.5, 4, 5, 6, 8, 10]
    cmap = nwsprecip()
    cmap.set_over('k')

    mp.pcolormesh(lons, lats, pcpn, clevs, cmap=cmap, latlon=True,
                  units='inch')
    lons, lats, vals, labels = get_data()
    mp.drawcounties()
    mp.plot_values(lons, lats, vals, "%s", labels=labels,
                   labelbuffer=1, labelcolor='white')

    mp.drawcities(labelbuffer=5, minarea=0.2)
    mp.postprocess(filename='test.png')

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def main():
    """Go Main"""
    grbs = pygrib.open('ds.snow.bin')
    # skip 1-off first field
    total = None
    lats = lons = None
    for grb in grbs[1:]:
        if lats is None:
            lats, lons = grb.latlons()
            total = grb['values']
            continue
        total += grb['values']
    # TODO tz-hack here
    analtime = grb.analDate - datetime.timedelta(hours=5)

    mp = MapPlot(
        sector='custom', west=-100, east=-92, north=45, south=41,
        axisbg='tan',
        title=("NWS Forecasted Accumulated Snowfall "
               "thru 7 PM 12 April 2019"),
        subtitle='NDFD Forecast Issued %s' % (
            analtime.strftime("%-I %p %-d %B %Y"), )
    )
    cmap = nwssnow()
    cmap.set_bad('tan')
    mp.pcolormesh(
        lons, lats, total * 39.3701,
        [0.01, 1, 2, 3, 4, 6, 8, 12, 18, 24, 30, 36],
        cmap=cmap,
        units='inch')

    mp.drawcounties()
    mp.drawcities()
    mp.postprocess(filename='test.png')
    mp.close()

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def test_drawcities():
    """Draw Cities"""
    mp = MapPlot(title='Fill and Draw Cities',
                 subtitle="This is my subtitle",
                 continentalcolor='blue',
                 sector='iowa', nocaption=True)
    mp.drawcities()
    return mp.fig

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def test_colorbar():
    """Run tests against the colorbar algorithm"""
    mp = MapPlot(sector='iowa', nocaption=True)
    cmap = plot.maue()
    cmap.set_under('white')
    clevs = list(range(0, 101, 10))
    norm = mpcolors.BoundaryNorm(clevs, cmap.N)
    mp.drawcities()
    mp.draw_colorbar(clevs, cmap, norm)
    return mp.fig

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def main():
    """Go Main"""
    nc = netCDF4.Dataset('/tmp/sfav2_CONUS_2018093012_to_2019021312.nc')
    lats = nc.variables['lat'][:]
    lons = nc.variables['lon'][:]
    data = nc.variables['Data'][:] * 1000. / 25.4
    nc.close()

    mp = MapPlot(
        sector='iowa', continentalcolor='tan',
        title=("National Snowfall Analysis - NOHRSC "
               "- Season Total Snowfall"),
        subtitle='Snowfall up until 7 AM 13 Feb 2019')
    cmap = plt.get_cmap('terrain_r')
    levs = [0.1, 2, 5, 8, 12, 18, 24, 30, 36, 42, 48]
    mp.pcolormesh(
        lons, lats, data, levs, cmap=cmap, units='inch', clip_on=False,
        spacing='proportional'
    )
    mp.drawcounties()
    mp.drawcities()
    mp.postprocess(filename='test.png')

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def main():
    """Go Main Go."""
    pgconn = get_dbconn('iem')
    cursor = pgconn.cursor()

    cursor.execute("""SELECT ST_x(geom) as lon, ST_y(geom) as lat,
    pday from summary_2016 s JOIN stations t on (s.iemid = t.iemid)
    where day = '2016-08-24' and network in ('WI_COOP', 'MN_COOP', 'IA_COOP')
    and pday > 0 ORDER by pday DESC""")
    llons = []
    llats = []
    vals = []
    for row in cursor:
        llons.append(row[0])
        llats.append(row[1])
        vals.append("%.2f" % (row[2], ))

    pgconn = get_dbconn('postgis')
    cursor = pgconn.cursor()

    cursor.execute("""SELECT ST_x(geom) as lon, ST_y(geom) as lat,
    max(magnitude) from lsrs_2016
    where wfo in ('DMX', 'DVN', 'ARX') and typetext = 'HEAVY RAIN' and
    valid > '2016-08-23' GROUP by lon, lat ORDER by max DESC""")
    for row in cursor:
        llons.append(row[0])
        llats.append(row[1])
        vals.append("%.2f" % (row[2], ))

    img = Image.open("p24h_201608241200.png")
    data = np.flipud(np.asarray(img))
    # 7000,3500 == -130,-60,55,25 ===  -100 to -90 then 38 to 45
    sample = data[1800:2501, 3000:4501]
    sample = np.where(sample == 255, 0, sample)
    data = sample * 0.01
    data = np.where(sample > 100, 1. + (sample - 100) * 0.05, data)
    data = np.where(sample > 180, 5. + (sample - 180) * 0.2, data)
    lons = np.arange(-100, -84.99, 0.01)
    lats = np.arange(38, 45.01, 0.01)

    x, y = np.meshgrid(lons, lats)

    buff = 0.5
    m = MapPlot(sector='custom', projection='aea', west=-93.2,
                east=-90.3, south=42.5,
                north=44.,
                title='NOAA MRMS 24 Hour RADAR-Only Precipitation Estimate',
                subtitle=("MRMS valid 7 AM 23 Aug 2016 to 7 AM 24 Aug 2016, "
                          "NWS Local Storm + COOP Reports Overlaid"))
    clevs = [0.01, 0.1, 0.25, 0.5, 0.75, 1, 1.5, 2, 2.5, 3, 3.5, 4, 5, 6, 8,
             10]

    m.contourf(x[:, :-1], y[:, :-1], data, clevs, cmap=nwsprecip())

    nt = NetworkTable("IA_ASOS")
    lo = []
    la = []
    va = []
    for sid in nt.sts.keys():
        lo.append(nt.sts[sid]['lon'])
        la.append(nt.sts[sid]['lat'])
        va.append(nt.sts[sid]['name'])

    # m.plot_values(lo, la, va, fmt='%s', textsize=10, color='black')
    m.drawcounties(zorder=4, linewidth=1.)
    m.drawcities(labelbuffer=25, textsize=10, color='white',
                 outlinecolor='#000000')
    m.textmask[:, :] = 0
    m.plot_values(llons, llats, vals, fmt='%s', labelbuffer=5)

    m.postprocess(filename='test.png')

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def main(argv):
    """Go Main Go"""
    v = argv[1]
    agg = argv[2]
    ts = datetime.date(2008, 1, 1)
    ts2 = datetime.date(2017, 12, 31)
    scenario = 0

    # suggested for runoff and precip
    if V2UNITS[v] in ['mm', 'inches']:
        colors = ['#ffffa6', '#9cf26d', '#76cc94', '#6399ba', '#5558a1']
    # suggested for detachment
    elif v in ['avg_loss', 'avg_loss_metric']:
        colors = ['#cbe3bb', '#c4ff4d', '#ffff4d', '#ffc44d', '#ff4d4d',
                  '#c34dee']
    # suggested for delivery
    elif v in ['avg_delivery', 'avg_delivery_metric']:
        colors = ['#ffffd2', '#ffff4d', '#ffe0a5', '#eeb74d', '#ba7c57',
                  '#96504d']
    cmap = mpcolors.ListedColormap(colors, 'james')
    cmap.set_under('white')
    cmap.set_over('black')

    pgconn = get_dbconn('idep')

    title = "for %s" % (ts.strftime("%-d %B %Y"),)
    if ts != ts2:
        title = "between %s and %s" % (ts.strftime("%-d %b %Y"),
                                       ts2.strftime("%-d %b %Y"))
    mp = MapPlot(axisbg='#EEEEEE', nologo=True, sector='iowa',
                 nocaption=True,
                 title=("DEP %s %s %s"
                        ) % (V2NAME[v.replace("_metric", "")],
                             "Yearly Average" if agg == 'avg' else 'Total',
                             title),
                 caption='Daily Erosion Project')

    df = read_postgis("""
    WITH data as (
      SELECT huc_12, extract(year from valid) as yr,
      sum(""" + v.replace("_metric", "") + """)  as d from results_by_huc12
      WHERE scenario = %s and valid >= %s and valid <= %s
      GROUP by huc_12, yr),

    agg as (
      SELECT huc_12, """ + agg + """(d) as d from data GROUP by huc_12)

    SELECT ST_Transform(simple_geom, 4326) as geo, coalesce(d.d, 0) as data
    from huc12 i LEFT JOIN agg d
    ON (i.huc_12 = d.huc_12) WHERE i.scenario = %s and i.states ~* 'IA'
    """, pgconn, params=(scenario, ts, ts2, scenario), geom_col='geo',
                      index_col=None)
    df['data'] = df['data'] * V2MULTI[v]
    if df['data'].max() < 0.01:
        bins = [0.01, 0.02, 0.03, 0.04, 0.05]
    else:
        bins = np.array(V2RAMP[v]) * (10. if agg == 'sum' else 1.)
    norm = mpcolors.BoundaryNorm(bins, cmap.N)

    # m.ax.add_geometries(df['geo'], ccrs.PlateCarree())
    for _, row in df.iterrows():
        c = cmap(norm([row['data'], ]))[0]
        arr = np.asarray(row['geo'].exterior)
        points = mp.ax.projection.transform_points(ccrs.Geodetic(),
                                                   arr[:, 0], arr[:, 1])
        p = Polygon(points[:, :2], fc=c, ec='k', zorder=2, lw=0.1)
        mp.ax.add_patch(p)

    mp.drawcounties()
    mp.drawcities()
    lbl = [round(_, 2) for _ in bins]
    u = "%s, Avg: %.2f" % (V2UNITS[v], df['data'].mean())
    mp.draw_colorbar(bins, cmap, norm,
                     clevlabels=lbl, units=u,
                     title="%s :: %s" % (V2NAME[v], V2UNITS[v]))
    plt.savefig('%s_%s_%s%s.png' % (ts.year, ts2.year, v,
                                    "_sum" if agg == 'sum' else ''))

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
import matplotlib.pyplot as plt
import datetime
import numpy as np

nc = netCDF4.Dataset("/mesonet/data/iemre/2016_mw_mrms_daily.nc", 'r')
idx = iemre.daily_offset(datetime.date(2016, 7, 3))
grid = np.zeros((len(nc.dimensions['lat']),
                 len(nc.dimensions['lon'])))
total = np.zeros((len(nc.dimensions['lat']),
                  len(nc.dimensions['lon'])))
for i, x in enumerate(range(idx, idx-60, -1)):
    total += nc.variables['p01d'][x, :, :]
    grid = np.where(np.logical_and(grid == 0,
                                   total > 25.4), i, grid)

m = MapPlot(sector='iowa', title='NOAA MRMS Q3: Number of Recent Days till Accumulating 1" of Precip',
  subtitle='valid 4 July 2016: based on per calendar day estimated preciptation, GaugeCorr and RadarOnly products')
lon = np.append(nc.variables['lon'][:], [-80.5])
print lon
lat = np.append(nc.variables['lat'][:], [49.])
print lat
x, y = np.meshgrid(lon, lat)
cmap = nwsprecip()
m.pcolormesh(x, y, grid,
             np.arange(0, 53, 7), cmap=cmap, units='days')
m.drawcounties()
m.drawcities()
m.postprocess(filename='test.png')
m.close()
nc.close()

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def run(basets, endts, view):
    """Generate this plot for the given basets"""

    df = read_sql("""SELECT state,
        max(magnitude) as val, ST_x(geom) as lon, ST_y(geom) as lat
        from lsrs WHERE type in ('S') and magnitude >= 0 and
        valid > %s and valid < %s GROUP by state, lon, lat
        """, POSTGIS, params=(basets, endts), index_col=None)
    df['used'] = False
    df['textplot'] = True
    df.sort_values(by='val', ascending=False, inplace=True)

    # Now, we need to add in zeros, lets say we are looking at a .25 degree box
    mybuffer = 0.75
    newrows = []
    for lat in np.arange(reference.MW_SOUTH, reference.MW_NORTH, mybuffer):
        for lon in np.arange(reference.MW_WEST, reference.MW_EAST, mybuffer):
            df2 = df[(df['lat'] >= lat) & (df['lat'] < (lat+mybuffer)) &
                     (df['lon'] >= lon) & (df['lon'] < (lon+mybuffer))]
            if len(df2.index) == 0:
                newrows.append(dict(lon=(lon+mybuffer/2.),
                                    lat=(lat+mybuffer/2.),
                                    val=0, used=True, textplot=False))
                continue
            maxval = df.at[df2.index[0], 'val']
            df.loc[df2[df2['val'] > (maxval * 0.5)].index, 'used'] = True
            df.loc[df2[df2['val'] < (maxval * 0.5)].index, 'textplot'] = False
    dfnew = pd.DataFrame(newrows)
    df = df.append(dfnew)
    cdf = df[df['used']]
    tdf = df[df['textplot']]

    rng = [0.01, 1, 2, 3, 4, 6, 8, 12, 18, 24, 30, 36]
    cmap = nwssnow()
    m = MapPlot(sector='iowa', axisbg='white',
                title="Local Storm Report Snowfall Total Analysis",
                subtitle=("Reports past 12 hours: %s"
                          "" % (endts.strftime("%d %b %Y %I:%M %p"), )))
    m.contourf(cdf['lon'].values, cdf['lat'].values, cdf['val'].values, rng,
               cmap=cmap)
    m.drawcounties()
    m.plot_values(tdf['lon'].values, tdf['lat'].values, tdf['val'].values,
                  fmt='%.1f')
    m.drawcities()
    pqstr = "plot c 000000000000 lsr_snowfall.png bogus png"
    m.postprocess(view=view, pqstr=pqstr)
    m.close()

    # slightly different title to help uniqueness
    m = MapPlot(sector='iowa', axisbg='white',
                title="Local Storm Report Snowfall Total Analysis",
                subtitle=("Reports valid over past 12 hours: %s"
                          "" % (endts.strftime("%d %b %Y %I:%M %p"), )))
    m.contourf(cdf['lon'].values, cdf['lat'].values, cdf['val'].values, rng,
               cmap=cmap)
    m.drawcounties()
    m.drawcities()
    pqstr = "plot c 000000000000 lsr_snowfall_nv.png bogus png"
    m.postprocess(view=view, pqstr=pqstr)
    m.close()

    m = MapPlot(sector='midwest', axisbg='white',
                title="Local Storm Report Snowfall Total Analysis",
                subtitle=("Reports past 12 hours: %s"
                          "" % (endts.strftime("%d %b %Y %I:%M %p"), )))
    m.contourf(cdf['lon'].values, cdf['lat'].values, cdf['val'].values, rng,
               cmap=cmap)
    m.drawcities()
    pqstr = "plot c 000000000000 mw_lsr_snowfall.png bogus png"
    m.postprocess(view=view, pqstr=pqstr)
    m.close()

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
m = MapPlot(sector='custom', projection='aea', west=-93.2,
            east=-90.3, south=42.5,
            north=44.,
            title='NOAA MRMS 24 Hour RADAR-Only Precipitation Estimate',
            subtitle=("MRMS valid 7 AM 23 Aug 2016 to 7 AM 24 Aug 2016, "
                      "NWS Local Storm + COOP Reports Overlaid"))
clevs = [0.01, 0.1, 0.25, 0.5, 0.75, 1, 1.5, 2, 2.5, 3, 3.5, 4, 5, 6, 8,
         10]

m.contourf(x[:, :-1], y[:, :-1], data, clevs, cmap=nwsprecip())

nt = NetworkTable("IA_ASOS")
lo = []
la = []
va = []
for sid in nt.sts.keys():
    lo.append(nt.sts[sid]['lon'])
    la.append(nt.sts[sid]['lat'])
    va.append(nt.sts[sid]['name'])

# m.plot_values(lo, la, va, fmt='%s', textsize=10, color='black')
m.map.drawcounties(zorder=4, linewidth=1.)
m.drawcities(labelbuffer=25, textsize=10, color='white',
             outlinecolor='#000000')
m.textmask[:, :] = 0
m.plot_values(llons, llats, vals, fmt='%s', labelbuffer=5)



m.postprocess(filename='test.png')

# 需要导入模块: from pyiem.plot import MapPlot [as 别名]
# 或者: from pyiem.plot.MapPlot import drawcities [as 别名]
def plotter(fdict):
    """ Go """
    import matplotlib
    matplotlib.use('agg')
    from pyiem.plot import MapPlot, nwsprecip
    ctx = util.get_autoplot_context(fdict, get_description())
    ptype = ctx['ptype']
    sdate = ctx['sdate']
    edate = ctx['edate']
    if sdate.year != edate.year:
        return 'Sorry, do not support multi-year plots yet!'
    days = (edate - sdate).days
    sector = ctx['sector']

    idx0 = iemre.daily_offset(sdate)
    idx1 = iemre.daily_offset(edate) + 1
    ncfn = "/mesonet/data/iemre/%s_mw_mrms_daily.nc" % (sdate.year, )
    if not os.path.isfile(ncfn):
        return "No MRMS data for that year, sorry."
    nc = netCDF4.Dataset(ncfn, 'r')
    lats = nc.variables['lat'][:]
    lons = nc.variables['lon'][:]
    if (idx1 - idx0) < 32:
        p01d = distance(np.sum(nc.variables['p01d'][idx0:idx1, :, :], 0),
                        'MM').value('IN')
    else:
        # Too much data can overwhelm this app, need to chunk it
        for i in range(idx0, idx1, 10):
            i2 = min([i+10, idx1])
            if idx0 == i:
                p01d = distance(np.sum(nc.variables['p01d'][i:i2, :, :], 0),
                                'MM').value('IN')
            else:
                p01d += distance(np.sum(nc.variables['p01d'][i:i2, :, :], 0),
                                 'MM').value('IN')
    nc.close()

    if sdate == edate:
        title = sdate.strftime("%-d %B %Y")
    else:
        title = "%s to %s (inclusive)" % (sdate.strftime("%-d %b %Y"),
                                          edate.strftime("%-d %b %Y"))
    if sector == 'midwest':
        state = None
    else:
        state = sector
        sector = 'state'
    m = MapPlot(sector=sector, state=state, axisbg='white', nocaption=True,
                title='NOAA MRMS Q3:: %s Total Precip' % (title,),
                subtitle='Data from NOAA MRMS Project, GaugeCorr and RadarOnly'
                )
    if np.ma.is_masked(np.max(p01d)):
        return 'Data Unavailable'
    clevs = [0.01, 0.1, 0.3, 0.5, 1, 1.5, 2, 2.5, 3, 3.5, 4, 5, 6, 8, 10]
    if days > 6:
        clevs = [0.01, 0.3, 0.5, 1, 1.5, 2, 3, 4, 5, 6, 7, 8, 10, 15, 20]
    if days > 29:
        clevs = [0.01, 0.5, 1, 2, 3, 4, 5, 6, 8, 10, 15, 20, 25, 30, 35]
    if days > 90:
        clevs = [0.01, 1, 2, 3, 4, 5, 6, 8, 10, 15, 20, 25, 30, 35, 40]
    x, y = np.meshgrid(lons, lats)
    cmap = nwsprecip()
    cmap.set_over('k')
    if ptype == 'c':
        m.contourf(x, y, p01d, clevs, cmap=cmap, label='inches')
    else:
        m.pcolormesh(x, y, p01d, clevs, cmap=cmap, label='inches')
    if sector != 'midwest':
        m.drawcounties()
        m.drawcities()

    return m.fig