Use case: anticipated profits based on precipitations probability¶

# Importation of Python modules 
import json
from datetime import datetime, timedelta
import re
import warnings

# The following modules must first be installed to use 
# this code out of Jupyter Notebook
import matplotlib.pyplot as plt
import matplotlib.dates as mdates
import numpy
from owslib.wms import WebMapService
import pandas
from tabulate import tabulate

# Ignore warnings from the OWSLib module
warnings.filterwarnings('ignore', module='owslib', category=UserWarning)

# Parameters choice
# Layer:
layer = 'REPS.DIAG.3_PRMM.ERGE5'
# Coordinates:
y, x = 49.288, -123.116
# Local time zone (in this exemple, the local time zone is UTC-07:00):
time_zone = -7

# bbox parameter
min_x, min_y, max_x, max_y = x - 0.25, y - 0.25, x + 0.25, y + 0.25

# WMS service connection
wms = WebMapService(
    'https://geo.weather.gc.ca/geomet?layers=REPS.DIAG.3_PRMM.ERGE5',
    version='1.3.0'
)

# Extraction of temporal information from metadata
def time_parameters(layer):
    start_time, end_time, interval = (wms[layer]
                                      .dimensions['time']['values'][0]
                                      .split('/')
                                      )
    iso_format = '%Y-%m-%dT%H:%M:%SZ'
    start_time = datetime.strptime(start_time, iso_format)
    end_time = datetime.strptime(end_time, iso_format)
    interval = int(re.sub(r'\D', '', interval))
    return start_time, end_time, interval


start_time, end_time, interval = time_parameters(layer)

# To use specific starting and ending time, remove the #
# from the next lines and replace the start_time and
# end_time with the desired values:
# start_time = 'YYYY-MM-DDThh:00'
# end_time = 'YYYY-MM-DDThh:00'
# fmt = '%Y-%m-%dT%H:%M'
# start_time = datetime.strptime(start_time, fmt) - timedelta(hours=time_zone)
# end_time = datetime.strptime(end_time, fmt) - timedelta(hours=time_zone)

# Calculation of date and time for available predictions
# (the time variable represents time at UTC±00:00)
time = [start_time]
local_time = [start_time + timedelta(hours=time_zone)]
while time[-1] < end_time:
    time.append(time[-1] + timedelta(hours=interval))
    local_time.append(time[-1] + timedelta(hours=time_zone))

# Loop to carry out the requests and extract the probabilities
def request(layer):
    info = []
    pixel_value = []
    for timestep in time:
        # WMS GetFeatureInfo query
        info.append(wms.getfeatureinfo(layers=[layer],
                                       srs='EPSG:4326',
                                       bbox=(min_x, min_y, max_x, max_y),
                                       size=(100, 100),
                                       format='image/jpeg',
                                       query_layers=[layer],
                                       info_format='application/json',
                                       xy=(50, 50),
                                       feature_count=1,
                                       time=str(timestep.isoformat()) + 'Z'
                                       ))
        # Probability extraction from the request's results
        text = info[-1].read().decode('utf-8')
        json_text = json.loads(text)
        pixel_value.append(
            float(json_text['features'][0]['properties']['value'])
            )
    
    return pixel_value

pixel_value = request(layer)

# Umbrellas sold per day in average when there is less than 30% 
# chance that there will be a minimum of 5 mm of precipitations
base = 3

# Profits per umbrella
umbrella_profit = 10.00

# Slope calculation data
# When the probability of precipitations is 30%...
x1 = 30
# ... 10 umbrellas are sold each hour
y1 = 10
# When the probability of precipitations is 100%...
x2 = 100
# ... 30 umbrellas are sold each hour
y2 = 30

# Slope calculation
slope = ((y2-y1)/(x2-x1))

# Open hours
opening = '09:00:00'
opening = datetime.strptime(opening, '%H:%M:%S')
closing = '21:00:00'
closing = datetime.strptime(closing, '%H:%M:%S')

# Prediction times that are within the open hours
open_hours = []
for timestep in local_time:
    if (timestep.time() > opening.time() 
            and timestep.time() < closing.time()):
        open_hours.append(timestep)

# Number of umbrellas sold each day independently of meteorological conditions
opening_interval = (opening + timedelta(hours=interval)).time()
umbrella = []
for timestep in local_time:
    new_day = timestep.time() < opening_interval
    if (umbrella == [] or new_day) and timestep in open_hours:
        umbrella.append(base)      
    else:
        umbrella.append(0)

# Number of umbrellas sold and anticipated profits
# depending on precipitations probability
cumulative_profit = []
for index, probability in enumerate(pixel_value):
    # Equation to calculate the number of umbrellas sold per hour
    eq = y1 + round((probability - 30) * slope)
    # Equation to calculate the number of umbrellas sold between 2 predictions
    eq2 = eq * interval 
    if local_time[index] in open_hours and probability > 30:
        if index == 0:
            umbrella[index] = umbrella[index] + eq2
        else: 
            umbrella[index] = umbrella[index] + umbrella[index - 1] + eq2
    elif index != 0:
        umbrella[index] = umbrella[index] + umbrella[index-1]                            
    cumulative_profit.append(umbrella[index] * umbrella_profit)