How to use the metpy.units.units function in MetPy

def test_convert_units(test_var):
    """Test in-place conversion of units."""
    test_var.metpy.convert_units('degC')

    # Check that variable metadata is updated
    assert units(test_var.attrs['units']) == units('degC')

    # Make sure we now get an array back with properly converted values
    assert test_var.metpy.unit_array.units == units.degC
    assert_almost_equal(test_var[0, 0, 0, 0], 18.44 * units.degC, 2)

def test_apparent_temperature_windchill():
    """Test that apparent temperature works when a windchill is calculated."""
    temperature = -5. * units.degC
    rel_humidity = 50. * units.percent
    wind = 35. * units('m/s')
    truth = -18.9357 * units.degC
    res = apparent_temperature(temperature, rel_humidity, wind)
    assert_almost_equal(res, truth, 0)

z = np.array([[5586387.00, 5584467.50, 5583147.50],
                  [5594407.00, 5592487.50, 5591307.50],
                  [5604707.50, 5603247.50, 5602527.50]]).T \
        * (9.80616 * units('m/s^2')) * 1e-3
    dx = np.deg2rad(0.25) * Re * np.cos(np.deg2rad(44))
    # Inverting dy since latitudes in array increase as you go up
    dy = -np.deg2rad(0.25) * Re
    f = (2 * omega * np.sin(np.deg2rad(44))).to('1/s')
    ug, vg = geostrophic_wind(z * units.m, f, dx, dy, dim_order='xy')
    true_u = np.array([[21.97512, 21.97512, 22.08005],
                       [31.89402, 32.69477, 33.73863],
                       [38.43922, 40.18805, 42.14609]])
    true_v = np.array([[-10.93621, -7.83859, -4.54839],
                       [-10.74533, -7.50152, -3.24262],
                       [-8.66612, -5.27816, -1.45282]])
    assert_almost_equal(ug[1, 1], true_u[1, 1] * units('m/s'), 2)
    assert_almost_equal(vg[1, 1], true_v[1, 1] * units('m/s'), 2)

# Convert data into a suitable format for metpy.
    _altitude = data[:,0] * units('m')
    p = mpcalc.height_to_pressure_std(_altitude)
    T = data[:,3] * units.degC
    Td = data[:,4] * units.degC
    wind_speed = data[:,1] * units('m/s')
    wind_direction = data[:,2] * units.degrees
    u, v = mpcalc.wind_components(wind_speed, wind_direction)


    fig = plt.figure(figsize=(6,8))
    skew = SkewT(fig=fig)
    skew.plot(p, T, 'r')
    skew.plot(p, Td, 'g')

    my_interval = np.arange(300, 1000, 50) * units('mbar')
    ix = mpcalc.resample_nn_1d(p, my_interval)
    skew.plot_barbs(p[ix], u[ix], v[ix])
    skew.ax.set_ylim(1000,300)
    skew.ax.set_xlim(-40, 30)
    skew.plot_dry_adiabats()

    heights = np.array([1, 2, 3, 4, 5, 6, 7, 8, 9]) * units.km
    std_pressures = mpcalc.height_to_pressure_std(heights)
    for height_tick, p_tick in zip(heights, std_pressures):
        trans, _, _ = skew.ax.get_yaxis_text1_transform(0)
        skew.ax.text(0.02, p_tick, '---{:~d}'.format(height_tick), transform=trans)

    plt.title("Sounding: " + title)

    if saveplot != None:
        fig.savefig(saveplot, bbox_inches='tight')

cmap = ctx["cmap"]
        if varname in ["rsds", "power_swdn"]:
            # Value is in W m**-2, we want MJ
            multi = (86400.0 / 1000000.0) if varname == "rsds" else 1
            data = nc.variables[varname][idx0, jslice, islice] * multi
            plot_units = "MJ d-1"
            clevs = np.arange(0, 37, 3.0)
            clevs[0] = 0.01
            clevstride = 1
        elif varname in ["wind_speed"]:
            data = (
                masked_array(
                    nc.variables[varname][idx0, jslice, islice],
                    units("meter / second"),
                )
                .to(units("mile / hour"))
                .m
            )
            plot_units = "mph"
            clevs = np.arange(0, 41, 2)
            clevs[0] = 0.01
            clevstride = 2
        elif varname in ["p01d", "p01d_12z", "snow_12z", "snowd_12z"]:
            # Value is in W m**-2, we want MJ
            data = (
                masked_array(
                    nc.variables[varname][idx0, jslice, islice], units("mm")
                )
                .to(units("inch"))
                .m
            )
            plot_units = "inch"

def update_precip(date, station, hdf):
    """Do the update work"""
    sts = utc(date.year, date.month, date.day, 7)
    ets = sts + datetime.timedelta(hours=24)
    ldf = hdf[
        (hdf["station"] == station)
        & (hdf["valid"] >= sts)
        & (hdf["valid"] < ets)
    ]

    newpday = ldf["precip_in"].sum()
    newpday_mm = (units("inch") * newpday).to(units("mm")).m
    set_iemacces(station, date, newpday)
    # update isusm
    pgconn = get_dbconn("isuag")
    cursor = pgconn.cursor()
    # daily
    cursor.execute(
        "UPDATE sm_daily SET rain_mm_tot_qc = %s, rain_mm_tot_f = 'E' "
        "WHERE valid = %s and station = %s",
        (newpday_mm, date, station),
    )
    for _, row in ldf.iterrows():
        # hourly
        LOG.debug(
            "set hourly %s %s %s", station, row["valid"], row["precip_in"]
        )
        total = (units("mm") * row["precip_in"]).to(units("inch")).m

(
                    "/mesonet/ARCHIVE/data/%Y/%m/%d/model/ffg/"
                    "5kmffg_%Y%m%d%H.grib2"
                )
            )
            if os.path.isfile(testfn):
                fn = testfn
                break
        if fn is None:
            raise NoDataFound("No valid grib data found!")
        grbs = pygrib.index(fn, "stepRange")
        grb = grbs.select(stepRange="0-%s" % (hour,))[0]
        lats, lons = grb.latlons()
        data = (
            masked_array(grb.values, data_units=units("mm"))
            .to(units("inch"))
            .m
        )
        plot.pcolormesh(lons, lats, data, bins, cmap=cmap)
        if ilabel:
            plot.drawcounties()
        df = pd.DataFrame()
    return plot.fig, df

""",
        pgconn,
        params=(station,),
        index_col=None,
    )
    if df.empty:
        raise NoDataFound("No Data was found.")
    # saturation vapor pressure
    # Convert sea level pressure to station pressure
    df["pressure"] = mcalc.add_height_to_pressure(
        df["slp"].values * units("millibars"),
        ctx["_nt"].sts[station]["elevation"] * units("m"),
    ).to(units("millibar"))
    # Compute the mixing ratio
    df["mixing_ratio"] = mcalc.mixing_ratio_from_relative_humidity(
        df["relh"].values * units("percent"),
        df["tmpf"].values * units("degF"),
        df["pressure"].values * units("millibars"),
    )
    # Compute the saturation mixing ratio
    df["saturation_mixingratio"] = mcalc.saturation_mixing_ratio(
        df["pressure"].values * units("millibars"),
        df["tmpf"].values * units("degF"),
    )
    df["vapor_pressure"] = mcalc.vapor_pressure(
        df["pressure"].values * units("millibars"),
        df["mixing_ratio"].values * units("kg/kg"),
    ).to(units("kPa"))
    df["saturation_vapor_pressure"] = mcalc.vapor_pressure(
        df["pressure"].values * units("millibars"),
        df["saturation_mixingratio"].values * units("kg/kg"),
    ).to(units("kPa"))

ds = ncss.get_data(hgt)

lon = ds.variables['lon'][:]
lat = ds.variables['lat'][:]

times = ds.variables[ds.variables['Geopotential_height_isobaric'].dimensions[0]]
vtime = num2date(times[:], units=times.units)


lev_500 = np.where(ds.variables['isobaric'][:] == 500)[0][0]

hght_500 = ds.variables['Geopotential_height_isobaric'][0, lev_500, :, :]
hght_500 = ndimage.gaussian_filter(hght_500, sigma=3, order=0) * units.meter

uwnd_500 = units('m/s') * ds.variables['u-component_of_wind_isobaric'][0, lev_500, :, :]
vwnd_500 = units('m/s') * ds.variables['v-component_of_wind_isobaric'][0, lev_500, :, :]

#######################################
# Begin Data Calculations
# -----------------------

dx, dy = mpcalc.lat_lon_grid_deltas(lon, lat)

f = mpcalc.coriolis_parameter(np.deg2rad(lat)).to(units('1/sec'))

avor = mpcalc.vorticity(uwnd_500, vwnd_500, dx, dy, dim_order='yx') + f

avor = ndimage.gaussian_filter(avor, sigma=3, order=0) * units('1/s')

vort_adv = mpcalc.advection(avor, [uwnd_500, vwnd_500], (dx, dy), dim_order='yx') * 1e9

#######################################

col_names = ['pressure', 'height', 'temperature', 'dewpoint', 'direction', 'speed']

df = pd.read_fwf(get_test_data('may4_sounding.txt', as_file_obj=False),
                 skiprows=5, usecols=[0, 1, 2, 3, 6, 7], names=col_names)

# Drop any rows with all NaN values for T, Td, winds
df = df.dropna(subset=('temperature', 'dewpoint', 'direction', 'speed'), how='all'
               ).reset_index(drop=True)

###########################################
# We will pull the data out of the example dataset into individual variables and
# assign units.

p = df['pressure'].values * units.hPa
T = df['temperature'].values * units.degC
Td = df['dewpoint'].values * units.degC
wind_speed = df['speed'].values * units.knots
wind_dir = df['direction'].values * units.degrees
u, v = mpcalc.wind_components(wind_speed, wind_dir)

###########################################
# Create a new figure. The dimensions here give a good aspect ratio.

fig = plt.figure(figsize=(9, 9))
add_metpy_logo(fig, 115, 100)
skew = SkewT(fig, rotation=45)

# Plot the data using normal plotting functions, in this case using
# log scaling in Y, as dictated by the typical meteorological plot.
skew.plot(p, T, 'r')
skew.plot(p, Td, 'g')