<jupyter_start><jupyter_text>2.2.c Temperature climatology in Spain (results)=================================================In this notebook we will finally generate the climatology results of daily temperature observations.Our final goal is still to know whether the temperature anomaly has influence on the power demand or not. Now, let's import the _classic stack_<jupyter_code>%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import warnings
warnings.filterwarnings('ignore')<jupyter_output><empty_output><jupyter_text>We do now need to load the compressed filed into memory.<jupyter_code>DATA = pd.read_csv("aemet_daily_observations_2009-2018.csv.bz2",compression='bz2')<jupyter_output><empty_output><jupyter_text>Let's have a quick look at the first and last lines of the file<jupyter_code>DATA.head(2)
DATA.tail(2)<jupyter_output><empty_output><jupyter_text>The observations cover the period end of 2009 end of 2018, 10 years of data. It seems enough to have a representative climatology.And data are already been thorugh some preparation.   Group by The **_groupby_** method is a very powerful tool for apply functions based on categories.We will use it to apply it to our _category_ ~ _day-of-the-year_.But, first of all, we add datetime properties to the date column and set it as index<jupyter_code>DATA.loc[:,'date'] = pd.DatetimeIndex(pd.to_datetime(DATA['date']))
DATA.set_index(['date'],inplace=True)
DATA.loc[:,'doy'] = DATA.index.strftime('Y%m%d')
DATA.head(2)
station_group = DATA.groupby(by='station')<jupyter_output><empty_output><jupyter_text>Let's loop throug the station_group object to have a first insight on what the groupby method has produce.<jupyter_code>for name,group_st in station_group:
    print ('%s (%s) %s' % (name,group_st.iloc[0].region ,len(group_st)))<jupyter_output>A CORUÑA (A CORUÑA) 3318
A CORUÑA AEROPUERTO (A CORUÑA) 3314
A POBRA DE TRIVES (OURENSE) 3194
ABLA (ALMERIA) 3257
ADRA (ALMERIA) 1672
ALAJAR (HUELVA) 3220
ALBACETE (ALBACETE) 3317
ALBACETE BASE AÉREA (ALBACETE) 3318
ALBORÁN (ALMERIA) 1711
ALCANTARILLA, BASE AÉREA (MURCIA) 3318
ALCAÑIZ (TERUEL) 3082
ALICANTE-ELCHE AEROPUERTO (ALICANTE) 3318
ALICANTE/ALACANT (ALICANTE) 3318
ALMERÍA AEROPUERTO (ALMERIA) 3318
ANDÚJAR (JAEN) 3087
ANTEQUERA (MALAGA) 2933
<NAME> (HUESCA) 2915
<NAME> (BURGOS) 3301
<NAME>   (NAVARRA) 2946
ARANJUEZ (MADRID) 3093
<NAME> (BARCELONA) 2975
<NAME> (JAEN) 3137
ASTURIAS AEROPUERTO (ASTURIAS) 3318
AUTILLA DEL PINO (PALENCIA) 3257
AYAMONTE (HUELVA) 3281
BADAJOZ AEROPUERTO (BADAJOZ) 3318
BARCELONA (BARCELONA) 3237
BARCELONA AEROPUERTO (BARCELONA) 3318
BARCELONA, FABRA (BARCELONA) 3318
BARDENAS REALES, BASE AÉREA  (NAVARRA) 3145
BAZA (GRANADA) 2890
BAZTAN, IRURITA  (NAVARRA) 2819
BELORADO (BURGOS) 2932
BENAVENTE (ZAMORA) [...]<jupyter_text>The next step is to groupby date of the year.<jupyter_code>for name,group_doy in group_st.groupby(by='doy'):
    print ('%s (%s) %4.2f %4.2f %4.2f' % (name,len(group_doy),np.mean(group_doy['tmin']),
                                          np.mean(group_doy['tavg']),np.mean(group_doy['tmax'])))<jupyter_output>Y0101 (9) 4.66 10.06 15.43
Y0102 (9) 4.17 10.10 16.04
Y0103 (9) 4.53 10.09 15.63
Y0104 (7) 4.84 10.77 16.69
Y0105 (7) 2.37 9.23 16.08
Y0106 (7) 3.20 9.37 15.56
Y0107 (7) 4.46 9.77 15.07
Y0108 (7) 3.57 9.54 15.50
Y0109 (7) 2.69 9.06 15.37
Y0110 (7) 4.27 10.19 16.07
Y0111 (7) 4.67 10.40 16.10
Y0112 (8) 4.00 9.74 15.50
Y0113 (8) 4.23 9.87 15.53
Y0114 (8) 3.39 9.01 14.67
Y0115 (8) 1.99 8.04 14.09
Y0116 (8) 3.34 8.50 13.67
Y0117 (9) 2.43 8.87 15.31
Y0118 (9) 3.53 8.91 14.31
Y0119 (7) 3.87 8.97 14.05
Y0120 (6) 1.62 7.48 13.33
Y0121 (9) 1.65 8.02 14.35
Y0122 (8) 3.07 9.01 14.97
Y0123 (9) 2.33 8.66 14.93
Y0124 (9) 2.85 9.46 16.05
Y0125 (9) 3.66 9.60 15.52
Y0126 (9) 3.77 9.91 16.11
Y0127 (9) 4.81 9.62 14.44
Y0128 (8) 3.46 8.94 14.40
Y0129 (8) 3.15 8.84 14.50
Y0130 (8) 3.24 9.41 15.60
Y0131 (9) 3.47 9.66 15.84
Y0201 (9) 2.79 9.52 16.23
Y0202 (9) 2.51 8.24 13.98
Y0203 (8) 2.45 8.55 14.62
Y0204 (8) 3.59 9.12 14.66
Y0205 (9) 3.54 9.52 15.50
Y0206 (9) 2.88 9.49 16.13
Y0207 (9) 3.10 9.14 15.20
Y0208 [...]<jupyter_text>Finally, all together<jupyter_code>line_of_result = []
for station,group_st in station_group:
    region = group_st.iloc[0].region
    for doy,group_doy in group_st.groupby(by='doy'):
        ndata = len(group_doy)
        line_of_result.append([station,region,ndata,doy,np.mean(group_doy['tmin']),
                                np.mean(group_doy['tavg']),np.mean(group_doy['tmax'])])
CLIMATOLOGY = pd.DataFrame.from_records(line_of_result)
CLIMATOLOGY.head()<jupyter_output><empty_output><jupyter_text>Let's add the correct columns<jupyter_code>CLIMATOLOGY.columns = ['station','region','ndata','doy','tmin','tavg','tmax']
CLIMATOLOGY.head()<jupyter_output><empty_output><jupyter_text>Interactive mode Let's have a quick interactive look<jupyter_code>import ipywidgets as widgets
from ipywidgets import interact, interact_manual,interactive, interactive_output
from ipywidgets import fixed, FloatSlider, Dropdown, HBox, Label, VBox, Layout
# Function to select region and station (based on the selected region...)

# Define the region widget (initialised at region = A CORUÑA)
region_list = sorted(list(set(CLIMATOLOGY['region'])))
region_widget = Dropdown(options=region_list,
                         value='A CORUÑA',
                         description='Region:'
                        )
                        
# Define the station widget (initialised at region = 'A CORUÑA' & statio = 'A CORUÑA')
station_widget = Dropdown(options=sorted(list(set(CLIMATOLOGY[CLIMATOLOGY['region']=='A CORUÑA']['station']))),
                         value='A CORUÑA',
                         description='Station:'
                         )
                        
# The upodate station list
def on_update_brand_widget(*args):
    station_widget.options = sorted(list(set(CLIMATOLOGY[CLIMATOLOGY['region']==region_widget.value]['station'])))

# The observe method to link station to region
region_widget.observe(on_update_brand_widget, 'value')

# Function
def plot_region_station(df,region,station,column):
    fig = plt.figure(figsize=(20,5))
    plt.grid()

    df_region = df.loc[df['region'] == region]
    DataFrame = df_region.loc[df_region['station'] == station]
    ndata = DataFrame['ndata'].values[0]
    plt.title('%s at %s (%s) ndata=%s' % (column,station,region,ndata),fontsize=10)
    display(plt.plot(DataFrame[column]))
Y = interactive(plot_region_station,df=fixed(CLIMATOLOGY),region=region_widget, station=station_widget,
               column=['tmin','tavg','tmax']);

controls_Left = VBox([Y.children[0],Y.children[1]])
controls_Right = VBox([Y.children[2]])

controls = HBox([controls_Left,controls_Right], layout = Layout(flex_flow='row wrap'))
output = Y.children[-1]
display(VBox([controls, output]))<jupyter_output><empty_output><jupyter_text>This climatology doesn't look like a soft climatology... 10 years may not be enough.In order to generate a typical climatology, we need to compute the mean over a window.<jupyter_code>window_widgets = widgets.IntSlider(value=5,min=1,max=24,step=2,description='Window:',disabled=False,
    continuous_update=False,orientation='horizontal',readout=True,readout_format='d')

doy_list = sorted(list(set(CLIMATOLOGY['doy'])))

# Function
def plot_window(df,region,station,column,longwindow):
    window = (longwindow-1)//2
    df_region = df.loc[df['region'] == region]
    df_station = df_region.loc[df_region['station'] == station]
    df_station.reset_index(inplace=True)
    line_of_result = []
    for doy in doy_list:
        doyindex = df_station[df_station['doy'] == doy].index[0]

        if (doyindex-window) < 0:
            df_concat = pd.concat([df_station.iloc[doyindex-window:],df_station.iloc[0:doyindex],
                                  df_station.iloc[doyindex:1+doyindex+window]])
        elif (doyindex+window) > 365:
            df_concat = pd.concat([df_station.iloc[doyindex-window:],
                                  df_station.iloc[:window-(365-doyindex)]])
        else:
            df_concat = pd.concat([df_station.iloc[doyindex-window:doyindex],
                                       df_station.iloc[doyindex:1+doyindex+window]])

        line_of_result.append([doy,np.mean(df_concat[column])])    

    DataFrame = pd.DataFrame.from_records(line_of_result)
    fig = plt.figure(figsize=(20,5))
    plt.grid()

    plt.title('%s at %s (%s) - window : %s days' % (column,station,region,longwindow),fontsize=10)
    display(plt.plot(DataFrame[1]))
Z = interactive(plot_window,df=fixed(CLIMATOLOGY),region=region_widget, station=station_widget,
               column=['tmin','tavg','tmax'], longwindow=window_widgets);

controls_Left = VBox([Z.children[0],Z.children[1]])
controls_Right = VBox([Z.children[2],Z.children[3]])

controls = HBox([controls_Left,controls_Right], layout = Layout(flex_flow='row wrap'))
output = Z.children[-1]
display(VBox([controls, output]))<jupyter_output><empty_output><jupyter_text>Climatology results Taking the mean value over a window produces a more representative climatology.Given a selected window days, let's recaculate the climatology for all stations.<jupyter_code># I've choosed a window covering 21 days
longwindow = 21
station_group = CLIMATOLOGY.groupby(by='station')
window = (longwindow-1)//2
line_of_result = []
for station,df_station in station_group:
    df_station.reset_index(inplace=True)
    region = df_station.iloc[0].region
    
    for doy in doy_list:
        try:
            doyindex = df_station[df_station['doy'] == doy].index[0]

            if (doyindex-window) < 0:
                df_concat = pd.concat([df_station.iloc[doyindex-window:],df_station.iloc[0:doyindex],
                                      df_station.iloc[doyindex:1+doyindex+window]])
            elif (doyindex+window) > 365:
                df_concat = pd.concat([df_station.iloc[doyindex-window:],
                                      df_station.iloc[:window-(365-doyindex)]])
            else:
                df_concat = pd.concat([df_station.iloc[doyindex-window:doyindex],
                                           df_station.iloc[doyindex:1+doyindex+window]])

            line_of_result.append([station,region,doy,np.mean(df_concat['tmin']),
                                    np.mean(df_concat['tavg']),np.mean(df_concat['tmax'])])
        except:
            break

CLIMATOLOGY_WINDOW = pd.DataFrame.from_records(line_of_result)
CLIMATOLOGY_WINDOW.columns = ['station','region','doy','tmin','tavg','tmax']
CLIMATOLOGY_WINDOW = pd.DataFrame.from_records(line_of_result)<jupyter_output><empty_output><jupyter_text>The daily temperature climatology is now ready!Let's export the result to a CSV file.<jupyter_code>CLIMATOLOGY_WINDOW.to_csv('aemet_daily_temperature_climatology_2009-2018_window-21.csv.bz2',
                          compression='bz2',float_format='%1.1f',index=False)<jupyter_output><empty_output>
StarcoderdataJupyterStructuredCleanDedup