A Python-only library for estimating global radiation

The problem

Rough estimates of global solar radiation are useful in many cases, like for example checking a radiometer’s response, or predicting evaporations-transpiration.

Computing these estimates is possible, using one of the many existing models. The pbl_met library, for example, allows doing this in two different ways (three, if we consider the old MPDA method).

Existing estimation software is most often written in Fortran (the pbl_met!), in C/C++ or other compiled language, which is fine if performance is at stake. There exist use cases, however, where performance is not that big issue, and maybe functionality is of greater interest; this is the area of use of Python and other interpreted languages. In them, support for numerical computations may be very strong, but specific physical modelling may be less available.

Library radest attempt filling this specific gap.

Using the radest library

Using radest is a simple task: you just have to import it. Wel’ll see soon how to do, using a tiny example.

For the library (and the example code) to work, your Python installation should support version 3.6 or more, and have a working installation of NumPy. For using the library on your own, an elementary knowledge of NumPy, in addition to plain Python, is recommended. In particular, an understanding of ‘datetime64’ and ‘timedelta64’, the NumPy equivalents of ‘datetime’ standard Python types, is highly helpful.

Nothing really difficult, by the way.

And here’s the example.

#!/usr/bin/env python3

import numpy as np
import radest
import os
import sys

if __name__ == "__main__":

    # Get parameters
    if len(sys.argv) != 9:
        print("radlist - Program generating a time series of global radiation estimates")
        print()
        print("Usage:")
        print()
        print("  ./radlist <Lat> <Lon> <Fuse> <Height> <Initial_ISO_DateTime> <Delta_Time> <Num_Steps> <Out_File>")
        print()
        print("where")
        print()
        print("  <Lat>, the latitude, is positive and increasing northwards (decimal degrees)")
        print("  <Lon>, longitude, is positive and increasing eastwards (decimal degrees)")
        print("  <Fuse> is an integer indicating the hours displacement from GMT")
        print("  <Height> is the height above ground (m)")
        print("  <Initial_ISO_DateTime> is a string like 2019-03-08T12:00:00, indicating a date and time")
        print("  <Delta_Time> is the (intgeer!) time step (s)")
        print("  <Num_Steps> is the number of time steps (greater than 0)")
        print("  <Out_File> is the output file name")
        print()
        print("Copyright 2019 by Mauri Favaron")
        print("This is open-source code, covered by the MIT license")
        print()
        sys.exit(1)
    lat   = float(sys.argv[1])
    lon   = float(sys.argv[2])
    fuse  = int(sys.argv[3])
    z     = float(sys.argv[4])
    first = np.datetime64(sys.argv[5], 's')
    delta = int(sys.argv[6])
    n     = int(sys.argv[7])
    out   = sys.argv[8]

    # Generate the time sequence corresponding to the desired time range
    tm = np.array([first + np.timedelta64(i*delta, 's') for i in range(n)])

    # Estimate global radiation
    ra = radest.ExtraterrestrialRadiation(tm,delta, lat, lon, fuse)
    rg = radest.GlobalRadiation(ra, z)

    # Print results
    f = open(out, "w")
    f.write("date, Rg\n")
    for i in range(len(tm)):
        f.write("%s, %f\n" % (str(tm[i]), rg[i]))
    f.close()

Some examples from Milan

Using radest allows to estimate the clear-sky global radiation. What would we then expect on a year, 2019 for example?

Here’s an answer:

Global radiation, estimated using radest at Milan position (Latitude 45.5°, Longitude 9.16°, Altitude 123m) on year 2019

The preceding plot shows quite well the envelope of maximal radiation one could expect on a day. But, what happens on a specific day? We can see in the next plot, where data from approximate season beginnings are shown.

Global radiation estimated by radest in Milan at four representative days: 21 June (red), 21 Mar and 21 Sep (pink) and 21 Dec (blue) 2019; notice the maxima do not occur exactly at 12:00:00, because of the combined effects of Milan longitude differing by 5.5° by the time zone central longitude, and the Time Equation

Is the time equation visible in radest?

The preceding plot suggests something is in fact discernible in radest data. But, what its magnitude is, and what’s its time evolution? In other term: if we measure the shift of the maximum radiation occurrence from the local “noon”, do we see something really similar to the Time Equation?

In order to check this, I’ve built a tiny application (which, incidentally, shows some data processing in NumPy/Python.

This is the code:

#!/usr/bin/env python3

import numpy as np
import radest
import os
import sys

if __name__ == "__main__":

    # Get parameters
    if len(sys.argv) != 6:
        print("timeq - Program generating an estimate of the Time Equation")
        print()
        print("Usage:")
        print()
        print("  ./timeq <Lat> <Lon> <Fuse> <Year> <Out_File>")
        print()
        print("where")
        print()
        print("  <Lat>, the latitude, is positive and increasing northwards (decimal degrees)")
        print("  <Lon>, longitude, is positive and increasing eastwards (decimal degrees)")
        print("  <Fuse> is an integer indicating the hours displacement from GMT")
        print("  <Year> is an integer designating the year for which processing should occur")
        print("  <Out_File> is the output file name")
        print()
        print("Copyright 2019 by Mauri Favaron")
        print("This is open-source code, covered by the MIT license")
        print()
        sys.exit(1)
    lat   = float(sys.argv[1])
    lon   = float(sys.argv[2])
    fuse  = int(sys.argv[3])
    year  = int(sys.argv[4])
    out   = sys.argv[5]

    # Main loop: Iterate over all days in dsired year, and process them
    initial_date = np.datetime64("%4.4d-01-01T00:00:00" % year, 's')
    final_date   = np.datetime64("%4.4d-01-01T00:00:00" % (year+1), 's')
    num_days     = ((final_date - initial_date) / np.timedelta64(1, 'D')).astype(int)
    days         = np.array([np.datetime64('2000-01-01', 'D') for i in range(num_days)])
    time_eq      = np.array([0.0 for i in range(num_days)])
    for day in range(num_days):

        start_of_day = initial_date + np.timedelta64(day, 'D')
        print(start_of_day)

        # Compute extraterrestrial radiation for all seconds in this day
        time_stamp = np.array([start_of_day + np.timedelta64(sec, 's') for sec in range(86400)])
        ra         = radest.ExtraterrestrialRadiation(time_stamp, 1, lat, lon, fuse)

        # Find the index corresponding to the maximum extraterrestrial radiation
        # and use it to calculate the time stamp at which it occurs
        idx_max_ra  = np.argmax(ra)
        time_max_ra = time_stamp[idx_max_ra]

        # Compute difference from noon, and save it
        noon = start_of_day + np.timedelta64(12, 'h')
        days[day]    = start_of_day
        time_eq[day] = (time_max_ra - noon) / np.timedelta64(1, 's')

    # Print results
    f = open(out, "w")
    f.write("date, tm.eqn\n")
    for i in range(len(days)):
        f.write("%s, %f\n" % (str(days[i]), time_eq[i]))
    f.close()

And this is the plot:

Time equation, as visible in Milan; the mean value, -1402 s, accounts for the difference between the time zone central meridian and site actual meridian. Changes over year take into account the Earth orbit eccentricity and the acceleration around Winter solstice.

As we can see, there is something similar to the Time Equation.

Other interesting effects

One of the thing the radiation estimate is able to do is, as with any mathematical model, eliciting questions, and providing an instrument for answering them – to the degree of accuracy permitted by the model.

And these two are questions, maybe useful in a high school sciences course:

  1. How does radiation change with latitude, as we moe instantaneously over a given meridian?
  2. Let’s go back to Milan: how does radiation change with altitude on a given day?

The first question can be solved using code which scans all longitudes on a given hour. The code is a bit trickier than usual, to cope with change of local hour on different fuses:

#!/usr/bin/env python3

import numpy as np
import radest
import os
import sys

if __name__ == "__main__":

    # Get parameters
    if len(sys.argv) != 1:
        print("radtst - Program generating various interesting tables for global radiation estimates")
        print()
        print("Usage:")
        print()
        print("  ./radtst")
        print()
        print("Copyright 2019 by Mauri Favaron")
        print("This is open-source code, covered by the MIT license")
        print()
        sys.exit(1)

    # Build data sets and print them, in sequence
    lat    = 0.
    lon    = 0.
    fuse   = 0
    z      = 0.
    time03 = np.array([np.datetime64('2019-03-21T00:00:00', 's')])
    time06 = np.array([np.datetime64('2019-06-21T00:00:00', 's')])
    time09 = np.array([np.datetime64('2019-09-21T00:00:00', 's')])
    time12 = np.array([np.datetime64('2019-12-21T00:00:00', 's')])
    rg03   = np.array([0. for ang in range(360)])
    rg06   = np.array([0. for ang in range(360)])
    rg09   = np.array([0. for ang in range(360)])
    rg12   = np.array([0. for ang in range(360)])
    out    = 'Equator_2019-03-21.csv'
    for ang in range(360):
        lon  = ang
        fuse = int(ang/24)
        if fuse == 360:
            fuse = 0
        local03 = time03 + np.timedelta64(fuse, 'h')
        local06 = time06 + np.timedelta64(fuse, 'h')
        local09 = time09 + np.timedelta64(fuse, 'h')
        local12 = time12 + np.timedelta64(fuse, 'h')
        ra  = radest.ExtraterrestrialRadiation(local03, 3600, lat, lon, fuse)
        rg03[ang] = radest.GlobalRadiation(ra, z)
        ra  = radest.ExtraterrestrialRadiation(local06, 3600, lat, lon, fuse)
        rg06[ang] = radest.GlobalRadiation(ra, z)
        ra  = radest.ExtraterrestrialRadiation(local09, 3600, lat, lon, fuse)
        rg09[ang] = radest.GlobalRadiation(ra, z)
        ra  = radest.ExtraterrestrialRadiation(local12, 3600, lat, lon, fuse)
        rg12[ang] = radest.GlobalRadiation(ra, z)

    # Print results
    f = open(out, "w")
    f.write("date, rg.03, rg.06, rg.09, rg.12\n")
    for i in range(len(rg03)):
        f.write("%d, %f, %f, %f, %f\n" % (i, rg03[i], rg06[i], rg09[i], rg12[i]))
    f.close()

Here is the plotted result:

Global radiation, estimated on June 21 (red), December 21 (blue), March 21 and September 21 (pink) at Equator, mean sea level, 00:00:00 UTC, and different longitudes. Note the maxima do not occur all at 180°: the shift, of seasonal nature, is due to the Equation of Time.

The second question is somewhat simpler to answer, as the fuse is now just the same. Here is the code:

#!/usr/bin/env python3

import numpy as np
import radest
import os
import sys

if __name__ == "__main__":

    # Get parameters
    if len(sys.argv) != 1:
        print("radtst2 - Program generating various interesting tables for global radiation estimates")
        print()
        print("Usage:")
        print()
        print("  ./radtst2.py")
        print()
        print("Copyright 2019 by Mauri Favaron")
        print("This is open-source code, covered by the MIT license")
        print()
        sys.exit(1)

    # Build data sets and print them, in sequence
    lat    = 0.
    lon    = 0.
    fuse   = 0
    z      = 0.
    time03 = np.array([np.datetime64('2019-03-21T12:00:00', 's')])
    time06 = np.array([np.datetime64('2019-06-21T12:00:00', 's')])
    time09 = np.array([np.datetime64('2019-09-21T12:00:00', 's')])
    time12 = np.array([np.datetime64('2019-12-21T12:00:00', 's')])
    lon = 0
    fuse = 0
    lats   = [i-89 for i in range(179)]
    rg03   = np.array([0. for ang in range(179)])
    rg06   = np.array([0. for ang in range(179)])
    rg09   = np.array([0. for ang in range(179)])
    rg12   = np.array([0. for ang in range(179)])
    out    = 'Equator_Latitudes.csv'
    for lat in lats:
        flat = float(lat)
        ra  = radest.ExtraterrestrialRadiation(time03, 3600, flat, lon, fuse)
        rg03[lat+89] = radest.GlobalRadiation(ra, z)
        ra  = radest.ExtraterrestrialRadiation(time06, 3600, flat, lon, fuse)
        rg06[lat+89] = radest.GlobalRadiation(ra, z)
        ra  = radest.ExtraterrestrialRadiation(time09, 3600, flat, lon, fuse)
        rg09[lat+89] = radest.GlobalRadiation(ra, z)
        ra  = radest.ExtraterrestrialRadiation(time12, 3600, flat, lon, fuse)
        rg12[lat+89] = radest.GlobalRadiation(ra, z)

    # Print results
    f = open(out, "w")
    f.write("lat, rg.03, rg.06, rg.09, rg.12\n")
    for i in range(len(rg03)):
        f.write("%d, %f, %f, %f, %f\n" % (i-89, rg03[i], rg06[i], rg09[i], rg12[i]))
    f.close()

And this is the corresponding plot:

Radiation change with latitude (longitude = 0°, height = 0m; red line for 21. 06. 2019, blue for 21. 12. 2019, pink lines for 21. 03. 2019 and 21. 09. 2019)

Conclusions

So, the radest library and its routines are here!

For you to use, if you desire.

The MIT license I’ve adopted is really permissive: you are free to use my software in both open-source and commercial applications, basically with no obligation to thank me.

But of course, a citation, or a public mention, of me and the library, is welcome and, from your side, a sign of kindness.

If you have questions, both related to the physics and use, you may ask to me, maybe even by posting a comment to this article.

Have a nice (and, of course, radiation-filled) day.

La radiazione solare globale massima attesa a Milano

“Radiazione solare globale”?

Una delle grandezze maggiormente utili, e di difficile misura, in ambiente urbano è la “radiazione solare globale”. Vediamo, intanto, cos’è e come la si misura.

Per definizione, la radiazione solare “globale” è la parte di spettro solare che arriva in un punto del suolo.

La maggior parte di questa radiazione è costituita da luce visibile, ma non mancano importanti contributi nelle parti infrarossa ed ultravioletta.

La radiazione solare globale si misura con il piranometro, di cui vediamo un esemplare nella figura qui sotto.

Esempio di radiometro globale: il modello SPP (Standard Precision Pyranometer) di Eppley Laboratory Inc.

Radiazioni “diretta” e “diffusa”

Per come sono costruite le sue cupole, il pirometro “conteggia” la radiazione globale che gli arriva da tutte le parti.

Secondo le indicazioni del Rapporto no.8 WMO, il piranometro andrebbe posizionato in modo da non risentire di ombreggiature e riflessioni – ma di questo problema delicato vedremo più avanti.

In primissima approssimazione, di giorno ed a cielo sgombro, le cose vanno circa così:

 

Ciò che vede un piranometro: condizioni ottimali: radiazione diretta e diffusa

La radiazione “riflessa”

Quando le condizioni non sono ideali, e ci sono nubi nel cielo, tutto dipende da quante sono queste, e messe dove.

Se la copertura nuvolosa è totale, o comunque nella direzione del Sole, allora  la radiazione diretta è quasi del tutto interrotta, ed al suolo giunge un’enorme quantità di radiazione diffusa.

Se ci sono poche nubi nel cielo, e poste in modo da non interrompere l’arrivo della radiazione diretta, assistiamo alla comparsa di una riflessione: dalle nubi, verso il radiometro.

Ciò che vede un piranometro: condizioni non ottimali per la presenza di nubi nel cielo; radiazione diretta però non intercettata; compare una nuova quota riflessa dalle nubi

In casi come questo, la radiazione che giunge al suolo (e, quindi, la lettura del piranometro) è più alta di quella che sperimenteremmo a cielo sereno. La radiazione in più si manifesta nei dati istantanei, dove di solito è visibile per “picchi”, di durata proporzionale al tempo che la nube ha impiegato a sfilare nel cielo.

Vrbanitas…

Nei siti urbani, tutto è più complicato.

A quanto abbiamo già visto si sommano ombreggiature e riflessioni da parte degli edifici e degli altri manufatti urbani, oltre che dagli alberi (sperando ve ne siano), animali…

Ciò che vede un piranometro: caso urbano

“Ma la notte, no!”

Per noi, abitanti dei Paesi temperati, la “notte” è quella parte del giorno durante la quale “il Sole non c’è”.

Il passaggio tra dì e notte, in realtà, avviene in modo graduale attraverso fasi di crepuscolo più o meno avvertibili.

Indubbiamente, comunque, a cavallo della mezzanotte in condizioni naturali assisteremmo al “profondo della notte”. E lì, decisamente, il Sole nel cielo non si vede.

Ciò vuol forse dire che non c’è luce, di notte? No:

Sulla superficie della Terra il buio assoluto è un fenomeno rarissimo, praticamente inesistente. Anche nella notte più oscura c’è un certo chiarore, dovuto alla presenza delle stelle nel cielo, ed alla fosforescenza rilasciata, alla superficie, da un’infinità di creature. In casi particolari, poi, ci sono fonti artificiali (lampade, incendi, …) e naturali (vulcani in eruzione, …)

La radiazione solare globale massima al suolo di Milano (se Milano non esistesse…)

Con ciò che abbiamo visto, l’idea di stimare la radiazione globale massima al suolo sembra di una complessità non comune.

E pure, esistono numerosi modelli adatti allo scopo, più o meno sofisticati.

Un esempio è il modello ASCE, reperibile in [ASCE_2005] in due forme, una indicativa ed una maggiormente accurata.

Sono, questi, modelli “non rifrattivi”, che cioè non tengono conto della rifrazione esercitata dalla distribuzione verticale dell’atmosfera.

La loro applicazione permette di ottenere una stima, che dipende dal tempo e dalla posizione (in misura minore dalla pressione, temperatura ed umidità relativa).

Si tratta, appunto, di una stima, e non di una previsione. Si fonda infatti sull’ipotesi, molto aggressiva, che il cielo sia completamente sgombro da nubi e che non vi siano riflessioni.

Vale, inoltre, al livello del mare.

Ipotizzando una temperatura (all’incirca annua!) di 20 °C, un’umidità relativa “annua” del 70%, ed una pressione di 1013 hPa, se calcoliamo la stima della radiazione globale e ne prendiamo i valori massimi giornalieri, otteniamo questo grafico (per il 2019):

Andamento dei massimi giornalieri della radiazione globale che vi sarebbe in assenza di nubi e di riflessioni. La curva superiore è stata calcolata in atmosfera limpida, mentre la inferiore con una concentrazione di polveri rappresentativa del Nord Italia

La natura di stima fa sì che le due curve del grafico siano puramente indicative. E comunque, vediamo che a Milano, a inizio Maggio, non fossimo in ambiente urbano, in una bella giornata potremmo aspettarci tra circa 780 e 820 W/mq (molto probabilmente più vicino a 780 che ad 820 a causa della presenza pervasiva di polvere nell’atmosfera.

Una vera tecnica di validazione?

I modelli di radiazione, che tentano di riprodurre un fenomeno molto complicato riducendolo alla sola posizione del Sole in un’atmosfera standardizzata priva di nubi, danno inevitabilmente solo indicazioni di prima approssimazione.

Le figure all’inizio dell’articolo, che mostrano l’andamento in condizioni reali, indicano alcuni altri elementi che andrebbero considerati, nell’ipotesi di un modello maggiormente sofisticato, od anche soltanto al fine di validare le misure.

Quello della validazione delle misure radiometriche, comunque, è un problema, e di non poco conto. Il piranometro, alla fin fine, è composto da una schiera di termocoppie, sensori notoriamente delicati, affacciate a superfici di differente colore o tonalità.

Le caratteristiche di questi toni, o colori, possono cambiare nel tempo, e così le caratteristiche delle termocoppie e della catena di acquisizione. Quest’ultima è particolarmente a rischio, i segnali elettrici prodotti dal piranometro essendo molto piccoli (dell’ordine di 15 microvolt per Watt al metro quadro).

Molto difficile, quindi, che il coefficiente di risposta si possa considerare costante nel tempo. Per conoscerne il valore, occorre una opportuna procedura di calibrazione.

Ma, rispetto a che cosa? Ad un sensore in campo? Magari, per quanto di riferimento, mai calibrato da anni? Oppure…

Insomma: una brutta gatta da pelare – per la quale gli istituti meteorologici ed i centri di ricerca senza dubbio hanno definito procedure opportune (alcuni produttori, ad esempio Eppley, offrono servizi di calibrazione).

Di fronte a questa difficoltà di ordine pratico, mi chiedo, ha senso validare i dati di radiazione? Oppure, è più percorribile la via di compiere dei controlli di plausibilità più o meno approfonditi, ma indicativi?

Riferimenti

[ASCE_2005] R.G. Allen, I.A. Walter, R.L. Elliott, T.A. Howell, D. Itenfisu, M.E. Jensen, R.L. Snyder (eds),  The ASCE Standard Reference Evapotranspiration Equation, ASCE, ISBN 9780784408056

 

 

Quanto sono davvero aperti gli “open data”?

Dati “aperti”: cosa sono?

Oggi si fa un gran parlare di open data: concetto molto importante e prezioso, in un mondo nel quale forte è la pressione a mettere i dati sotto tutela, e soprattutto quelli sulla cui base si possono prendere decisioni importanti per la salute o il reddito dei cittadini.

Ma cosa sono, di preciso, i “dati aperti”, almeno qui in Italia? Una risposta, pragmatica, è: quelli resi disponibili dalla Pubblica Amministrazione e da altri soggetti, coperti dalla licenza IODL (Italian Open Data License), il cui testo aggiornato è disponibile alla URL https://www.dati.gov.it/content/italian-open-data-license-v20 .

Ad una prima lettura (superficiale) della licenza, mi sembra di poter dire che il suo scopo principale è di rendere i dati liberamente accessibili da parte di chiunque.

In questo articolo cercheremo di comprendere se il libero accesso è, anche, accompagnato da libero uso, e in caso negativo quali ne siano gli ostacoli, e le possibili contromisure.

Open-data come open-source?

Esistono altre licenze utilizzate per garantire la diffusione libera di dati (ad esempio la Creative Commons), me tutte si rifanno, a quanto comprendo, ad un sottofondo comune, simile, a prima vista, a quello per il software open source.

Simile, ma non identico. C’è infatti una differenza importantissima tra un software disponibile in forma di sorgente, e “dati” numerici.

Nel primo caso, il sorgente costituisce una sorta di specifica di ultima istanza, non ambigua e deterministicamente associabile al software in questione. Si può supporre che almeno una parte degli utenti finali abbia le competenze, il tempo e la voglia per comprendere ill funzionamento della cosa, giudicarne della sua adeguatezza ai propri casi in base alla propria esperienza, ed eventualmente apportare le correzioni del caso.

Con i dati, nulla di tutto ciò è possibile. C’è un’asimmetria ineliminabile di potere, in questo caso, tra chi quei dati li ha prodotti, e chi li usa. Con chi, poi, li ha raccolti e resi disponibili nel mezzo.

Se non si fa attenzione a ciò che si fa, la diffusione di open data potrebbe, così, costituire uno dei più colossali ed efficaci veicoli di propagazione di fake news.

Un esempio, tra i tantissimi esistenti

Come esempio di dati aperti, prendiamo le temperature rilevate in città nel 2018 rese disponibili dal Comune di Milano alla URL http://dati.comune.milano.it/dataset/ds686-rilevazione-temperature-anno-2018.

Come spesso accade la pagina di accesso non si limita ai link per lo scarico dei dati, ma contiene una tabella di metadati. Questi contribuiscono a rendere fruibili i dati, o dovrebbero farlo.

Dalla pagina dei metadati vediamo che:

  • La copertura dei dati va dall’1 Dicembre 2019 al 12 Dicembre 2019 (dato che oggi siamo al 7 Aprile 2019, presumo che queste due date siano nate da un refuso, e che i dati vadano dal primo di Gennaio 2018 al 31 Dicembre 2018.
  • L’editore del “data set” è il Comune di Milano.
  • Il titolare dei dati è ARPA Lombardia.
  • Il creatore dei dati è il SIAD – Unità Sviluppo Open Data, sempre del Comune di Milano a giudicare dalla partita IVA.

Più in basso, nella pagina dei metadati, vediamo che la loro origine è tacciabile alla pagina di scarico del sito di ARPA Lombardia, di cui è riportata la URL.

Quanto ai dati veri e propri, selezionando l’opzione “download” nel bottone di scarico possiamo farli comparire in una pagina testuale, dalla quale possiamo selezionarli, copiarli nella schermata di un editor di testo, e salvarli in un archivio CSV, per usarli con calma.

Queste le prime righe:

Zone;Id Sensore;Data-Ora;Media
Lambrate;2001;01/01/2018 00:00;3,5
Zavattari;5920;01/01/2018 00:00;3,2
Juvara;5909;01/01/2018 00:00;4
Feltre;8162;01/01/2018 00:00;3,2
Brera;5897;01/01/2018 00:00;3,6
Marche;5911;01/01/2018 00:00;3,8
Lambrate;2001;01/01/2018 01:00;3,3
Zavattari;5920;01/01/2018 01:00;3,1
Juvara;5909;01/01/2018 01:00;3,8
Feltre;8162;01/01/2018 01:00;3,2
Brera;5897;01/01/2018 01:00;3,4
Marche;5911;01/01/2018 01:00;3,5

Ci sono tutti. Ed è confortante vedere che, in effetti, partono dalle ore 00 del giorno 01. 01. 2019, e non dal 01. 12. 2019 come dichiarato nella pagina dei metadati.

Il fatto che i dati siano tutti presenti, però, porta con sé anche una difficoltà: per il modo in cui sono ordinati nel file, ogni riga una singola stazione, non è possibile utilizzarli direttamente, ad esempio per farne una serie storica.

E le righe sono tante. In un anno ci sono 8760 ore (24 in più, se bisestile). Per sei stazioni, il totale arriva a 52560.  Troppe, per poterle prendere in esame una ad una.

Sarebbe decisamente meglio, se invece che in un singolo file contenente tutti i dati di tutte le stazioni, vi fosse invece un file per stazione.

Ma se desideriamo porre i dati in questa nuova, e più utilizzabile forma, dobbiamo fare un certo lavoro. Per esempio, scrivere e usare uno script R come questo:

get.comune <- function() {

# Read Milan's open data to a data frame, and complement # it with time stamps
d <- read.csv2("Temp_2018_ComuneMilano.csv",
   stringsAsFactors = FALSE)
yr <- substring(d$Data.Ora, first=7, last=10)
mo <- substring(d$Data.Ora, first=4, last=5)
dy <- substring(d$Data.Ora, first=1, last=2)
tm <- paste(substring(d$Data.Ora, first=12),
   "00", sep=":")
date.time <- paste(yr,"-",mo,"-",dy," ",tm, sep="")
d$Time.Stamp <- as.POSIXct(date.time, tz="UTC")

# Remove now-unused columns, and arrange the 
# remaining ones in better user-order
d <- d[c(5,1,2,4)]
names(d) <- c("Time.Stamp", "Zone", "Sensor.ID", "Temp")

# Partition data to (for example) separate files,
# one per station
stations <- unique(d$Zone)
for(station.name in stations) {
  file.name <- paste("Temp_", station.name, "_2018.csv", sep="")
  e <- d[d$Zone==station.name,]
  e$Zone <- NULL
  write.csv2(e, file=file.name, row.names=FALSE)
}

# Leave (yielding back the transformed data set,
# for debug if needed)
return(d)

}

Piccoli (ma forse grandissimi) problemi concreti

Interessante, dover fare del lavoro non banale, per accedere a dei dati “aperti”. Il cittadino che desiderasse accedere ai dati di via Zavattari dovrebbe assumere un perito informatico, oppure imparare a programmare in R.

Ma in effetti nella pratica professionale accade sempre di dover “pulire” dati di provenienza esterna, prima di riuscire ad utilizzarli. Un professionista, però, ha dalla sua strumenti di lavoro, conoscenze, pratica: cioè, il risultato di lunghi ed ingenti investimenti. Nel caso, il concetto di dati “aperti” risulta un ossimoro: se al professionista i dati risultano non soddisfacenti, avrà sempre modo se necessario di “scassinarli” e, con opportune martellate informatiche, rimetterli in squadra.

L’utente finale, però, questa possibilità non ce l’ha.

Forse, l’origine del problema va ricercata in una mancanza di comunicazione e, da parte di chi rilascia le informazioni e detiene grazie a ciò una posizione di potere, assenza di empatia.

Intanto, tra necessità proprie e mentalità di chi i dati li ha preparati (e che magari non è nemmeno una persona umana) intercorre un abisso insondabile di non-detto, ma solo, nell’ipotesi migliore, immaginato, in modo quasi certamente divergente tra le due parti.

La responsabilità di colmare questo abisso sta a chi sa di più: la parte debole, nell’interazione, non sa di non sapere, e non ha alcuno strumento concettuale per porsi domande, concettuali od operative.

E i dati del Comune di Milano non fanno eccezione.

Colpevole, dunque, il Comune di Milano? Di sicuro, di non aver pensato abbastanza, prima di diffondere informazioni che potrebbero essere usate contro di lui.

Specchio, però, di un malvezzo molto diffuso, radicato in una profonda ignoranza dei mezzi informatici, a sua volta tacciabile all’arretratezza cultural del nostro Paese. Basti pensare alla valanga di siti, di Banche o della Pubblica Amministrazione, realizzati in modo incomprensibile. Come la pagina di una nota banca che, per indurre i Clienti a passare ad una nuova tecnologia di autenticazione, permette dalla stessa schermata di scaricare un’app per Android, o per iOS, o per Windows, pretendendo che gli utenti (spesso persone anziane sole, sappiano cos’hanno tra le mani, e che l’applicazione va scaricata sul telefono, invece che sul computer: questo, ai tempi di HTML5 e dei siti “adattivi”).

Il problema tutto sommato veniale incontrato per caso sugli open data del Comune di Milano ne solleva forse un altro, di portata immensamente maggiore e più pericolosa: di una classe dirigente che usa a piene mani le tecnologie informatiche, senza comprenderne minimamene pericoli ed opportunità, e senza alcuna possibilità di svolgere il suo compito di guida. Come, se ricordate, decenni e decenni fa, i responsabili del Fronte Polisario (i Tuareg!) che avevano acquistato dei caccia Mirage, salvo constatare che sulla sabbia del deserto non potevano atterrare. O che,  casomai lo facessero, restavano impantanati e non ripartivano più. Con il risultato di una spesa ingente, che non ha comportato alcun risultato apprezzabile: quegli aerei non hanno fatto alcun danno al nemico, e se per quello non hanno mai nemmeno consegnato una cartolina postale.

Il Principio della Trota.

Sono le ore 15:21 del 28 Aprile 1979, e, evento importantissimo nella nostra storia, il signor M.R. ha pescato un’alborella. Una sola, per quanto alla foce del Tresa tra Luino e Germignaga ai tempi qualcuna ce ne fosse ancora.

Alle 16:30, l’alborella nel cestino era ancora da sola, e contava tristissima le alghe, che la Tresa di tanto in tanto portava verso il Lago Maggiore.

Alle 17, il signor M.R. decise che ne aveva abbastanza di dover spiegare come mai il cestino fosse ancora così desolatamente casi vuoto, liberò la sua preda prima che morisse di tristezza, e si incamminò verso il bar.

E lì, racconta che ti racconta, l’alborella subì una misteriosa, quanto comune, trasformazione. Già nel racconto di M.R. era divenuta due volte più grande che nella realtà. E, forzutissime.

Ma il suo amico S.D. (un po’ sordo) capì “trota bella” invece che “alborella”, con il che il pesciolino raddoppiò ancora la sua lunghezza.

Al quarto passaggio, l’alborella era divenuta un mostro marino, un capodoglio forse, che per un misterioso accidente della Natura, e sospinto da un’ammirevole forza di volontà, aveva superato (con una bella shakerata e qualche danno per l’orgoglio) la turbina Kaplan dell’impianto di Tornavento, risalendo il Ticino sino a Somma, Angera, ed oltre…

Questo fenomeno, che in società è molto comune, e che nell’arte della pesca assume una curiosa proclività verso l’aumentare e mai diminuire, si chiama “distorsione semantica”, ed avviene ogni volta che tra le due parti impegnate nel trasferimento di informazioni intercorre una differenza esperienziale. Maggiore la distanza, più evidenti le distorsioni.

Se poi andiamo ad analizzare in dettaglio cosa accade (e proprio tra i dettagli, com’è noto, si nasconde il Diavolo), vediamo che qualche volta le informazioni passano, ma ripetute con le parole di chi le ha ricevute, diverse nel senso da quelle di chi le ha emesse. Qualche altra volta, informazioni che chi riceve considera irrilevanti, e che magari non lo erano, vengono pietosamente omesse. Altre volte ancora, chi riceve ritiene l’informazione ricevuta non adatta, pur avendone compreso il senso, e a ristruttura in una forma che ritiene più comoda. Oppure ancora, chi riceve ha tra le mani qualcosa di incompleto, ma non ha modo di sapere dell’incompletezza, e diffonde i dati ricevuti come fossero completi ed esaustivi.

Per dirla in breve: persino informazioni “oggettive” come dati misurati corrono il rischio d’essere sottoposti ad un pesante, ed inconsapevole, processo di invenzione creativa che ne modifica il significato.

Adesso, che abbiamo finalmente il nostro script di ordinamento, possiamo finalmente dividere i dati per stazione, e verificare se qualche distorsione c’è stata anche nel caso dei nostri dati.

La stazione più vicina alla casa dei miei è quella di Milano Zavattari: decido di puntare l’attenzione su di essa, e come prima operazione ne conto i dati.

Sorpresa!

Le righe non sono 8760! Dovrebbero essere proprio tante così, dal momento che il 2018 non era bisestile, per cui 365 per 24 fa 8760 e…

Ma i dati che ho ricevuto sono di meno: 8737. Ne mancano ventitrè.

Una rapida scorsa ad uno dei file rivela il problema: i dati arrivano, in effetti, sino al 31 Dicembre 2018, ma, si fermano alle ore 00! Quelli dalle 01 alle 23 non ci sono proprio.

Chi ha predisposto i dati, con ogni probabilità, lo ha fatto usando la pagina di scarico di ARPA Lombardia, passando come date estreme

  2018-01-01,  2018-12-31

dimenticandosi che 2018-12-31 vuol dire in realtà 2018-12-31 00:00:00. Un errore classico, da neofiti, o da chi ha poco tempo. Gli estremi corretti erano (forse in modo poco intuitivo)

2018-01-01,  2019-01-01

Ma, ormai, tant’è: i dati “open” sono capitozzati irreversibilmente, almeno sino a quando il Comune di Milano provvederà alla sperata correzione.

Tra i metadati dell’apposita pagina compaiono, anche, i nomi della persona che ha materialmente predisposto il file con i dati, e del dirigente cui questa fa capo.

Fossimo gente di corte vedute ci verrebbe spontaneo gettar la croce addosso a queste due, addossando loro ogni colpa per questi, e presumibilmente altri, errori, e pretendendo, se non la correzione degli stessi, almeno immediata soddisfazione.

Ritengo però più costruttivo, ed interessante, affrontare il problema in modo più sistematico, avendo ravvisato in esso caratteri che trascendono il Comune di Milano, e gettano una luce (sinistra, in parte) sulla natura stessa degli “open data”.

Sono, questi “open data” del Comune di Milano, in effetti dati riportati. Che qualcuno ha preparato, e qualcun altro ha preso in carico, trasformato in base a fini che nei metadati non mi risultano dichiarati (magari mi sbaglio), impacchettato, e distribuito in modo più o meno conveniente secondo a propria visione soggettiva.

In quanto “open data”, il focus, l’attenzione, è andata sull’apertura, sulla ampia diffusione.

Però, all’apparenza, ha trascurato molti altri aspetti. Ha assunto una visione, per così dire, “a trapano”, lasciando perdere il contesto generale.

In questo modo, ha anche involontariamente indirizzato le inevitabili omissioni e distorsioni, che si verificano quando un insieme di dati sia trasformato e ri-distribuito, in una direzione di massima efficienza, ma, possibilmente, minimo senso – e, difficile fruibilità: vedi il mio script R.

Se volessimo andare alla radice del problema, dovremmo forse dirigere la nostra attenzione agli organi politici che hanno deciso di produrre gli “open data”, allocato risorse (insufficienti?), deliberato azioni, senza dare un preciso mandato tecnico a chi poi avrebbe dovuto fare, ed anzi, delegando a questi ogni scelta operativa, con la solita possibilità di facile sconfessione in caso di problemi. E, malvezzo comune della politica italiana, evitando accuratamente di lasciarsi indicare per nome e cognome nella pagina dei metadati.

Ma. Siamo sicur* di aver definito interamente i termini del problema?

Forse, no. Potrebbe esserci dietro qualcosa di molto più profondo, e meno controllabile, che una manciata di refusi facili da sanare.

Qualcosa che, in perfetta buona fede, ma in modo profondissimamente sbagliato, stiamo facendo tutti insieme, con le nostre azioni, ed anche, concedetemi, con le nostre aspettative. Qualcosa, che cozza contro una realtà di fatto, dura come al solito, radicata nelle leggi della fisica e della matematica.

Andiamo all’origine?

I dati “open” di temperatura del Comune di Milano vengono, per ammissione dello stesso, da ARPA Lombardia. Ed anche lì, in fondo, sono dati “open”: aperti, accessibilissimi, tramite una pagina di scarico di facile impiego (salvo il significato della data “finale”).

Prendiamo, allora, i dati della stazione ARPA di Zavattari, facendo l’interrogazione così come la condurrebbe una persona comune, e cioè chiedendo quelli orari, in formato CSV, dal giorno 2018-01-01 al 2018-12-31.

Dopo qualche minuto dalla richiesta, all’indirizzo di posta elettronica che abbiamo indicato giunge un messaggio contenente, in allegato, i dati richiesti.

Che sono, nel caso della temperatura, composti da due file: i dati veri e propri, e la legenda.

Vediamola, perché è interessante:

Legenda
-999 Valore mancante o invalido

Riepilogo estrazione
Id_Stazione,Nome_Stazione,Id_Sensore,Nome_Sensore,Unita_Misura,Periodo_Dal,Periodo_Al,Operatore
503,Milano P.zza Zavattari,5920,Temperatura,°C,2018/01/01 00:00,2018/12/31 00:00,Media

La seconda riga, interessantissima, ci dice che esiste un valore speciale che rappresenta i dati invalidi.

Questa informazione, non rilasciata esplicitamente dal Comune di Milano, che per i suoi “open data” ha scelto una via diversa: i dati invalidi sono lasciati vuoti.

Una piccola, insignificante differenza. Che, però,  ha reso i dati in teoria aperti, la cui origine è attribuita ad ARPA Lombardia, diversi dagli originali.

Le righe successive, invece, hanno una funzione per lo più contabile: la richiesta di dati avrebbe potuto riguardare più sensori, e qui avremmo trovato i riferimenti a ciascuno di essi, insieme al suo significato complessivo. Noi di sensore ne abbiamo richiesto uno solo, e la riga che gli compete ci può sembrare poco importante.

Veniamo ai dati. Quelli di ARPA Lombardia sono fatti così:

Id Sensore,Data-Ora, Media
5920,2018/01/01 00:00,3.2
5920,2018/01/01 01:00,3.1
5920,2018/01/01 02:00,2.8
5920,2018/01/01 03:00,2.4
5920,2018/01/01 04:00,2.5
5920,2018/01/01 05:00,2.8

Interessante. Qui, 2018/01/01. Nei dati “aperti” del Comune di Milano, 01/01/2018. La stessa data, espressa però in due forme diverse. Prova, questa, che il Comune di Milano non ha soltanto acquisito i dati ARPA per riversarli come sono. Li ha, invece, letti, magari filtrati, e riscritti in una forma certamente diversa. Ci sono altre piccole differenze: il separatore decimale (la virgola, all’italiana, invece che il punto; il punto e virgola come separatore di campo invece della virgola.

Quanto pesante è stato il filtraggio? Non ne ho idea. Nei metadati non se ne parla.

Se andiamo a prendere i primi valori generati dal mio script per Zavattari, troviamo queste prime righe:

"Time.Stamp";"Sensor.ID";"Temp"
2018-01-01 00:00:00;5920;3,2
2018-01-01 01:00:00;5920;3,1
2018-01-01 02:00:00;5920;2,8
2018-01-01 03:00:00;5920;2,4
2018-01-01 04:00:00;5920;2,5
2018-01-01 05:00:00;5920;2,8

I valori di temperatura coincidono, quindi, “ci siamo”.

Se guardiamo all’intera serie dei dati, vediamo che i valori dei dati validi coincidono con quelli ARPA. Il che, però, non dice alcunché di universale: un anno diverso, o un’altra stazione, e magari qualche differenza si vedrebbe. Sarebbe molto più comodo, però, se qualcuno si fosse preso la briga di confermare questa uguaglianza.

Se, poi, andiamo in fondo al file ARPA, troviamo queste righe:

5920,2018/12/30 19:00,12.1
5920,2018/12/30 20:00,11.2
5920,2018/12/30 21:00,11.0
5920,2018/12/30 22:00,10.6
5920,2018/12/30 23:00,8.1
5920,2018/12/31 00:00,6.1

E’ così confermato, che richiedendo i dati al sito ARPA sino al 31 Dicembre, l’insieme che otteniamo a quell’ora si ferma.

Il significato delle marche temporali

La legenda ci conferma su una cosa che avevamo intuito, ma che sino ad ora non avevo resa esplicita. Le temperature che compaiono nei file sono medie orarie.

Medie, di che? Di misure elementari. Quali, di preciso?

Se i numeri sono medie orarie, allora si riferiscono  ad un intervallo temporale lungo un’ora. Intervallo, a rigore, contraddistinto da due marche temporali, una di inizio e una di fine.

L’unica marca temporale che correda ogni dato di provenienza ARPA Lombardia si riferisce per convenzione all’istante finale dell’intervallo di mediazione.

A priori, avrebbe potuto essere qualsiasi altro valore entro l’ora: l’inizio, la metà, o, appunto, la fine. O, qualsiasi altro valore, possibilmente stabilito da una qualche regola.

Quindi, la marca 2018-01-01 00:00:00 si riferisce all’intervallo di tempo che va dalle ore 23:00:00 del giorno 2017-12-31 alle 00:00:00 del 2018-01-01. E così, per tutte le altre righe dell’insieme di dati.

Comode da generare (per chi costruisce i sistemi di acquisizione dei dati), le marche temporali posticipate sono però in grado di suscitare una certa confusione. Molte persone considerano più intuitive le marche temporali anticipate. Alcune altre (ed io tra di loro) preferirebbero come marca temporale il centro dell’intervallo di mediazione.

Tutte scelte legittime, con propri vantaggi e svantaggi. Non esiste alcuna ragione, gusto estetico soggettivo a parte, per preferire l’una all’altra. Ciò che mi preme osservare, è che non esiste consenso circa le marche temporali di periodo, cosicché, se non fosse disponibile l’informazione a priori che la marca dei dati è posticipata, si rischierebbe di compiere un errore grave nell’interpretazione della posizione dei dati nel tempo.

Se al momento della distribuzione dei dati “aperti” questa informazione chiave è omessa, la loro comprensione finisce minata, e l’uso più difficile.

Certo, un’utente molto capace, e molto fortunata, potrebbe notare che il massimo del giorno-tipo delle temperature non occorre a mezzodì, ma sfasato di un’ora, ed avrebbe modo di attribuire ai dati una marca temporale per lei più sensata.

Ma dovrebbe, appunto, essere capace e fortunata.

O disporre di informazioni da insider trader, come quella che vi ho data io.

Dietro un numero

Ambiguità del linguaggio naturale

I dati di temperatura, lo abbiamo già visto, sono “medie”.

Basta, questa parola, per caratterizzarne in modo esaustivo natura e proprietà dei nostri “open data”, e darci elementi che ci permettano a priori di prevederne il valore?

La risposta breve è, “No”.

Ma, vediamo i particolari.

Partiamo, intanto, da un’ipotesi di lavoro. Quando la legenda dei dati ARPA parla di “media”, immaginiamoci di intendere “media aritmetica”. Qualcosa, insomma, di simile a questa formula:

\[ \overline{x}=\frac{1}{n}\sum_{i=1}^{n}x_{i} \]

Questa è la definizione più comune del termine “media”, ed è ciò che ci immaginiamo quando pensiamo ai polli di Trilussa, e ricordiamo sin dalle scuole medie.

A priori ci sarebbero anche altri tipi di “media”. Ad esempio, la “media armonica”,

\[ x_{h}=\frac{n}{\sum_{i=1}^{n}\frac{1}{x_{i}} } \]

Oppure, la “media geometrica”:

\[ x_{g}=\left(\prod_{i=1}^{n}x_{i}\right)^\frac{1}{n} \]

Oppure, altri tipi di “media” più esoterici, come ad esempio l’uscita di un filtro autoregressivo nel dominio del tempo.

La mia ipotesi che per “media” debba proprio intendersi quella aritmetica affonda le sue radici nell’ambiguità intrinseca del linguaggio naturale, che, a differenza dei linguaggi artificiali più diffusi (ad esempio i linguaggi di programmazione) richiede in modo esplicito l’intervento di un’intelligenza che nella miriade di possibili significati sappia scegliere quello, o quelli, più sensati in base al caso.

A volte, come nel caso che stiamo analizzando, il contesto tecnico non da alcun appiglio ad una intelligenza puramente neutrale. Infatti, a priori non ci sarebbe alcuna ragione di preferire una media armonica, o geometrica, a quella aritmetica.

La mia scelta a favore della media aritmetica, in fondo, non si basa su dati puramente tecnici, ma sulla psicologia e sulla tradizione: da un lato, i meteorologi sono abituati alla media aritmetica, e storcerebbero il naso di fronte a qualcos’altro; dall’altro, le “medie” sono prodotte direttamente dai sistemi automatici di registrazione dati che compongono, insieme agli strumenti, una stazione meteorologica: e chi costruisce sistemi automatici di registrazione dati, quasi sempre adotta dandola per scontata e senza nemmeno dirne la “media aritmetica”.

A priori, non è detto che l’utente dei dati aperti” condivida questo background, e potrebbe così immaginarsi qualcosa d’altro. Aspettativa legittima, sintantoché nessuno gli dice in modo esplicito che media è stata realmente utilizzata.

“Ma che importa,” chiederà a questo punto qualcuno, “sapere com’è stata calcolata la media?”

Vero, rispondo. In alcuni casi, sapere come sono state ottenute le medie non è così importante. Precisamente, nei casi in cui il dato “non ha un’importanza critica”, ma solo un interesse, per così dire, indicativo e qualitativo.

Esistono, però, altri casi, nei quali il dato ha un’importanza critica. Ad esempio, nelle cause legali, nei contratti, nella ricerca scientifica, nella meccanica di precisione ed in altri innumerevoli campi dell’ingegneria. Ed altrove ancora.

Lì, sapere se i dati sono stati ottenuti usando una media aritmetica (quindi uno stimatore del valore atteso di un processo casuale che sta dietro alla manifestazione visibile delle “misure” che hanno formato i dati) piuttosto che di altro tipo (con qualche altra proprietà profonda desiderabile) non è per nulla indifferente: conoscere o meno questo dettaglio permette, o no, di fare sulle “medie” delle affermazioni “forti”, del tipo, “Le cose stanno così, e ti spiego perché.”

Oltre i dati “aperti”?

Quanto gli “open data” del Comune di Milano hanno evidenziato potrebbe indicare un fatto importante, ed una potenzialità preoccupante.

La diffusione di dati “aperti” è un’intrapresa essenzialmente informatica.

Utile, ma, ad oggi, non per tutti. Come abbiamo veduto, il loro utilizzo reale comporta l’impiego di mezzi sofisticati, e conoscenze non banali perché tutto ciò dia luogo ad interpretazioni sensate.

Cosa manca?

Secondo me, due cose importanti:

  • Metadati accurati ed esaustivi, e,
  • Cultura.

A queste, mi permetto di aggiungerne una terza, non meno importante, ma di carattere sistemico, e a differenza delle due qui sopra, di difficile soluzione:

  • Consapevolezza tecnica ed eccellenza da parte della classe dirigente.

Nel caso esempio specifico, ai metadati, e ad una formazione più accurata dei dati da distribuire, può senz’altro provvedere, con le risorse disponibili e con poco sforzo, il Comune di Milano.

Più in generale: ma siamo sicuri* che i dati “aperti” siano sufficienti?

Soprattutto adesso, con la valanga di numeri che l’Internet delle Cose potrebbe riversare su di noi. Sull’apertura, ho pochi dubbi che ci sarà.

Ma che dire, della qualità? Dell’usabilità?

Già oggi è difficile validare i dati acquisiti per via convenzionale: fondi e risorse scarseggiano. Ma domani, o dopodomani, sarà materialmente ancora possibile farlo?

Si parla, di questi tempi, di big data. Di “identificazione proattiva dei guasti”. Di machine learning.

Ma chi garantirà che i numeri che ci inonderanno abbiano un rapporto preciso e raccontabile della realtà?

La risposta, forse, è nell’andare oltre.

Dai dati “aperti”, ai dati trasparenti.

Corredati da metadati che dicano tutto ciò che c’è da sapere al dettaglio necessario, e che non si possano disgiungere dai dati.

Magari, ad un “livello di qualità” dei dati, con una classificazione del tipo “usabili da utenti finali”, “usabili con cautela da professionisti”, “potenziale spazzatura”.

Un’eventuale licenza sui “transparent data” afferirebbe, inevitabilmente, non solo ai dati in quanto tali, ma all’intera catena che li produce. Hardware e software inclusi.

A Servizi Territorio srl stiamo, ad esempio, lavorando ad una nuova famiglia di sensori intelligenti con diffusione dei dati tramite IoT. Rispetto ai sensori intelligenti comunemente reperibili sul mercato, questi pongono l’attenzione non sulla diffusione veloce di numeri, ma sulla distribuzione di dati di qualità definita.

Che si possano usare, ad esempio, in uno studio sull’isola di calore. O più in generale, quando i dati abbiano davvero un valore, e la natura della loro formazione non sia indifferente.

Molto ancora, però, c’è da fare. Incluso, cambiare la percezione comune, ed in certo modo elevarla rispetto allo stato attuale, sufficiente forse nel secolo scorso, ma non più adatto per sopravvivere in tempi globali ed incerti.

Cominciando, magari, a riflettere tutti (cittadin*-utenti inclus*) su una licenza per dati trasparenti, che siano davvero utilizzabili senza problemi da parte di (potenzialmente) chiunque.

 

La pbl_met cresce…

Un problema piccolo, ma concreto

Molte grandezze meteorologiche e micro-meteorologiche sono misurate in un numero finito di punti, ma definite in uno spazio a più dimensioni: la superficie della Terra, oppure l’atmosfera, oppure ancora gli interstizi del suolo.

E’ lì, e non solo nei punti di misura, che vorremmo conoscere lo stato delle cose. A chiederlo sono la logica elementare, ma anche (in casi specifici come quello dell’inquinamento dell’aria) la normativa.

Nasce così il problema, annoso, di “spazializzare le misure”, locuzione non proprio elegantissima, ma di grande valore pratico.

Un modo sensato (tra i tanti, sensati, possibili) per estendere un insieme di misure finito ad una regione dello spazio è di usare tutte le conoscenze disponibili, racchiuderle in un modello fisico, ed applicarlo alle misure – direttamente, o dopo una fase di assimilazione.

Nel caso, il modello fisico produce delle previsioni spaziali, che sono valutate in corrispondenza di una griglia di punti recettori nello spazio. Griglia ancora finita, ma talmente fitta da permettere di ricostruire una buona immagine della grandezza che ci interessa.

Il passo successivo è, appunto, ricostruire l’immagine.

Nella pratica, ciò può essere fatto campionando i risultati del modello fisico su una griglia ancora più fitta, usando questa volta un metodo di interpolazione (molto più economico del modello fisico, in termini computazionali). E su quel campione, applicare finalmente una delle tante tecniche di rendering delle immagini, come ad esempio la costruzione delle curve di iso-livello.

Nella pbl_met

Proprio in questo momento, la pbl_met si sta arricchendo di una nuova classe, TwoDimensionalField, deputata proprio al trattamento di campi bidimensionali, ricostruiti a partire da un insieme di punti sperimentali (che, incidentalmente, non sono costretti a stare su di una grigliaa regolare).

Il lavoro è in progresso, e farò sapere quando sarà finalmente completato.

Posso però dire, che i metodi che prevedo di implementare sono molto semplici (per ora):

  • Inverso del quadrato (o di una potenza sensata qualsiasi) della distanza.
  • Metodo di Cressmann.

In futuro, vedremo. In letteratura ci sono dei bellissimi metodi molto più promettenti di questi, ad esempio quello descritto in F. Uboldi, C. Lussana, M. Salvati, “Three-dimensional spatial interpolation of surface observations from high-resolution local networks”, Meteorol. Appl., 2008.

Vi terrò informat*. 🙂

Diario di una stazione (micro-)meteorologica urbana – Seconda puntata: Ma in città le cosa vanno davvero in modo del tutto diverso che fuori?

Molte ragioni di dubbio

La stazione sperimentale micro-meteorologica di Cinisello Balsamo è posta, come ho detto nella prima puntata, proprio in mezzo all’urban canopy. E non di un villaggio qualsiasi: no, di una città densissima, estremamente rugosa.

Come se non bastasse, il punto di misura è posto praticamente in centro città. Come possiamo vedere qui:

Coincidenza vuole, che nelle vicinanze della stazione micro-meteorologica sperimentale vi sia la stazione micro-meteorologica di Parco Nord Cinisello, parte della rete SHAKEUP, e di proprietà ARPA Lombardia.

Dico nelle vicinanze, ma le due stazioni operano comunque in condizioni molto diverse. Intanto, per la distanza: poco più di un chilometro e mezzo.

E poi, soprattutto, per il contesto nel quale le due stazioni sono collocate, che possiamo vedere dalla prossima immagine:

Guardando l’ultima figura, ed immaginandoci d’essere particelle d’aria passate dalla stazione di Parco Nord, per arrivare alla stazione sperimentale di Cinisello, dovremmo riconoscere d’essere piuttosto sballottate nel viaggio, passando da un contesto che assomiglia molto da vicino a quello rurale, ad uno pesantemente urbano. Vien voglia di credere che, qualunque informazione contenessimo al passaggio da Parco Nord, questa verrebbe perduta una volta raggiunta la stazione di Cinisello.

E’ davvero così?

Correlazioni sorprendenti

Il grafico che segue confronta la direzione del vento rilevata a Parco Nord con quella di Cinisello:

Confronto tra le direzioni del vento rilevate dalla stazione SHAKEUP di Cinisello Parco Nord e dalla stazione sperimentale di Cinisello Balsamo – Dati: Anno 2016

Ciò che vediamo non è una nube di punti più o meno rotonda, indicativa di una bassa o nulla correlazione: anzi. E’ assolutamente evidente una fortissima correlazione. Con una certa non-linearità nelle direzioni oltre i 150° da Nord.

Cosa vuol dire, questa scoperta?

Una cosa importante: il flusso d’aria sopra la città di Cinisello Balsamo, dentro l’urban canopy, reca ancora un’importante traccia della circolazione “non cittadina”. Un’impronta.

Che non sia identità, lo vediamo dal prossimo grafico:

Confronto tra le velocità del vento rilevate dalla stazione SHAKEUP di Cinisello Parco Nord e dalla stazione sperimentale di Cinisello Balsamo – Dati: Anno 2016

La maggior frizione dovuta alla forte rugosità della città ha frenato il vento, che vediamo aver perso più della metà di tutta la sua velocità “rurale” (e qualche volta quasi tutta).

Vien da chiedersi, se una situazione del genere (direzione fortemente correlata, velocità molto cambiata) si possa applicare a tutte le stazioni urbane vicine ad un nodo SHAKEUP. La mia risposta è… non lo so!

Immagino, questo sì, che le cose vadano diversamente città per città, sito di misura per sito di misura. E che quanto vediamo alla stazione sperimentale di Cinisello Balsamo sia la conseguenza d’una coincidenza fortunata.

Ma soprattutto, possiamo dire una cosa: il mondo reale è sempre pieno di sorprese. Di storie, che racconta con voce sommessa. Non cercare d’imporgli delle regole, della prevedibilità, più di quelle che hanno certamente senso: nel caso finiresti col parlare tu. E rischieresti di dimenticarti di ascoltare.

 

Dalle misure alle mappe. Prima puntata: Dati “spaziali”?!

Il problema

La forma dei dati

Sarà accaduto moltissime volte anche a voi: avete ricevuto una lista di misure, prese in punti diversi più o meno nello stesso istante, e vi hanno chiesto di estenderle ad un’intera area più o meno grande.

Benvenuti nel club: anche a voi è stato posto il problema di “spazializzare” i dati.

Dato che la cosa è accaduta anche a me, che mi ha dato alcuni problemi, e che soffrire senza un motivo valido mi sembra cosa tutto sommato poco astuta, desidererei condividere alcune esperienze e idee, così come le ho maturate in questi ultimi più-di-vent’anni.

In questa prima puntata vedremo insieme il problema, ed alcune sue sfaccettature. Nelle puntate che seguiranno vedremo, invece, alcuni metodi pratici, con relativi vantaggi e svantaggi.

Sembra facile, ma…

Mi ricordo, quando avevo dieci anni, un Carosello in cui un buffo personaggio diceva qualcosa del tipo “Si si si si. Sembra facile.” E so per certo che non si riferiva ad un modello di caffettiera Moa di un noto produttore. No: era un messaggio trasversale, che bucava il tempo, e si occupava della spazializzazione dei dati.

“Spazializzare” i dati sembra, a prima vista, un problema semplice. Tipo, l’analogo in tre dimensioni del tracciare a mano la linea continua che raccorda i punti sperimentali.

A complicare il problema, però, c’è un elemento: i punti sperimentali sono affetti da errore sperimentale. Quelli che in prima battuta ci si presentano (o ci vengono presentati) come numeri, in realtà a stretto rigore andrebbero immaginati alla stregua di variabili casuali, con una ben definita distribuzione e, possiamo stante certi, una dispersione intorno alla media maggiore di zero.

Non basta: c’è anche il problema della scelta dei punti, che molto spesso non è per nulla libera.

Un Tipico Problema Malposto

Appunto: malposto. La spazializzazione di un insieme di misure è un procedimento intrinsecamente “visivo”: viene bene, se il risultato soddisfa il nostro senso estetico. E, avessimo a disposizione del tempo, ed un piccolo laboratorio di scultura, non sarebbe difficile plasmare una bella superficie che interpola od approssima il nostro insieme di misure.

Non appena però tentiamo di fare la stessa cosa (rappresentare i dati tramite una superficie “bella”) con tecniche algoritmiche, intanto per cominciare otterremo un risultato che dipende dall’algoritmo usato.

E nessuna delle soluzioni risulta davvero soddisfacente.

Il che porta inevitabilmente a scegliere una soluzione di compromesso. Cioè, una forma funzionale, un modello matematico, che non dia risultati troppo brutti e, nello stesso tempo, che coaguli attorno a sé un minimo di consenso.

Per capire cosa intendo dire, prendiamo un caso semplice:

Grafico e mappa di “isoconcentrazione” di un “inquinante” la cui concentrazione alla superficie del suolo è descritta esattamente dall’equazione c=exp(-(x^2+y^2)/100), valutata su una griglia “fitta” di 129 x 129 punti: la situazione ideale

Niente di che. Una situazione di questo tipo si verifica spontaneamente quando (caso piuttosto comune) un inquinante diffonde in matrice fluida in un mezzo poroso uniforme, se il fluido se ne sta fermo. Cosa che accade, per esempio, in certe falde acquifere. Casi analoghi, meno frequenti, si incontrano anche in atmosfera.

Il grafico precedente è stato costruito immaginando di poter campionare la concentrazione a proprio piacere, senza alcuna considerazione per tempi e denaro, compiendo nel caso più di 16000 misure, tutte nello stesso preciso istante.

Però, tempo e denaro esistono, come esistono recinti, muri, canali, e in generale tutte le tracce lasciate dalla nostra specie su, e dentro, il territorio. Ciò porta a campionare i punti su una grigia fatta più o meno così:

Griglia irregolare di punti di campionamento della concentrazione

Usando lo stesso metodo di interpolazione adoperato per le più-di-sedicimila-misure, ma limitandoci alle 16 che vediamo qui sopra, la figura che otteniamo è decisamente diversa:

Grafico e mappa di “isoconcentrazione” di un “inquinante” la cui concentrazione alla superficie del suolo è descritta esattamente dall’equazione c=exp(-(x^2+y^2)/100), valutata su una griglia “irregolare” costituita da 16 punti: la situazione più comune

Lasciatemelo dire: questo ultimo grafico è molto più brutto del primo.

Indubbiamente, non abbiamo fatto un grande sforzo. Abbiamo applicato in modo brutale il metodo di interpolazione “inverso del quadrato della distanza” (il default del pacchetto software che ho usato per comporre le mappe), senza domandarmi nulla di eventuali errori di misura od altro.

E come spesso accade in questi casi, ho evocato l’effetto GIGOGarbage In (in questo caso un metodo di interpolazione scelto a casaccio), Garbage Out.

Sorge spontanea, allora, la domanda: è possibile fare di meglio?

La risposta breve è: Sì, si può. Impegnandocisi, e conoscendo i propri dati.

Per la risposta lunga, dovrete avere la pazienza di aspettare le prossime puntate…

 

Diario di una stazione (micro-)meteorologica urbana – Prima puntata: Nascita

Antefatto

Correva il… Mah. 2015?

Certe date la memoria fa fatica a ricordarle: non sono come le nascite vere. O, forse, sì? Nel caso, però, per così dire “in miniatura”.

Perché a Servizi Territorio avevamo bisogno di un “muletto”, una stazione di misura esattamente uguale a quelle che avremmo di lì a poco consegnato per un grosso lavoro, in modo da essere assolutamente sicuri che il sistema di acquisizione non contenesse errori 🙂    E soprattutto, una volta scoperti i “bachi” mancanti all’appello, avremmo avuto modo di correggerli.

Missione compiuta. E poi…

La stazione è rimasta lì, sul tetto dell’ufficio, con il suo bellissimo anemometro ultrasonico tri-assiale USA-1 della Metek GmbH, ed un esemplare “arricchito” del nostro sistema eddy covariance MeteoFlux Core V2 (scusate le parolacce)..

Posizione della stazione micro-meteorologica sperimentale di Cinisello Balsamo. Vista di dettaglio del punto di misura.
Una vista più in grande della collocazione della stazione micro-meteorologica sperimentale di Cinisello Balsamo

Con vergogna, devo confessare che me la sono dimenticata lì, la stazione…

La pressione del collaudo era finita, tutto bene, si poteva passare ad altre cose più urgenti, più o meno frenetiche.

Lei intanto, la stazione, acquisiva, elaborava, registrava…

Anni dopo, chiedendo al mio sistema di fare un elenco dei nodi di rete, ne noto uno il cui indirizzo mi diceva qualcosa…

Attimi di rossore. La stazione sul tetto!

Era rimasta là, abbandonata, lontana dal cuore, ma, accesa (consumando così poco che non c’eravamo accorti della sua esistenza in vita).

Tento di collegarmi.

Non ci posso credere.

Dati, dat, dati.

Una montagna.

Per tutto questo tempo, il mio cucciolo MeteoFlux Core ha compiuto indefessamente la sua missione. Povera stella, senza che nessuno glie lo chiedesse.

“Il bello dei sistemi embedded, quando funzionano,” mi dissi, “fanno le loro cose stando nell’ombra, senza farsi notare.”

E, col regalo di qualche gigabyte di dati grezzi sonici…

Un’occasione inattesa

Dire che il tetto dell’ufficio di Servizi Territorio sia un sito di misura ideale, ecco…

Però, siamo proprio in mezzo all’urban canopy, come dicono gli esperti veri: quella parte di atmosfera tra tetti, balconi, terrazze, vie. Sopra la città, ma non tanto da lasciar libera l’aria di fluire come meglio crede.

La canopea urbana (?) (spero la parola esista davvero, comunque, è “quella roba lì”) è, anche, sede di fenomeni importantissimi. E’ lì, che la città e l’aria scambiano energia, schifezze, e quant’altro. Ed è lì, che la turbolenza ed il flusso dell’aria assumono caratteristiche tutte loro.

Insomma: cercando di collaudare il nostro MeteoFlux Core, eravamo riusciti a costruire, senza averlo pianificato, una stazione micro-meteorologica urbana.

Micro-meteorologica?

Una stazione insolita

Le stazioni meteorologiche sul tetto, o sul balcone, abbondano, e dunque una più, una meno non fa poi una gran notizia.

Questa, però, è micro-meteorologica: non è proprio la stessa cosa.

Al centro, l’anemometro ultrasonico tri-assiale USA-1 della stazione micro-meteorologica urbana di Servizi Territorio srl. Più in basso a sinistra vediamo invece gli strumenti della stazione “Cinisello Balsamo” della rete urbana della Fondazione OMD (Osservatorio Meteorologico di Milano Duomo).
Da solo, l’anemometro ultrasonico USA-1 misura, ad altissima risoluzione, le tre componenti del vento e la temperatura “virtuale” (comprensiva degli effetti dell’umidità) dell’aria.

Vediamo, allora, di capirci qualcosa.

Cenni minimi di micro-meteorologia.

Il prefisso “micro” ci suggerisce che la micro-meteorologia non ha molto a vedere con le previsioni del tempo, od il moto delle grandi masse d’aria. Evoca, piuttosto, un concetto di “piccolo”. Ma: piccolo quanto?

Meglio che ci rifacciamo ad una definizione assestata.

Secondo Glickman (Glickman, 2000)

“La micro-meteorologia è la parte di meteorologia che si occupa delle osservazioni e dei processi che si svolgono alle scale più piccole dello spazio e del tempo, approssimativamente più piccole di 1 chilometro ed 1 giorno. I processi micro-meteorologici si limitano a sottili strati che risentono dell’influenza diretta della frizione (fenomeni a scala leggermente maggiore come le termiche convettive non sono parte della micro-meteorologia). Quindi, il soggetto della micro-meteorologia è la parte più bassa dello Strato Limite Planetario (PBL), cioè, in particolare, lo Strato Superficiale. I processi di scambio di energia, gas, eccetera, tra l’atmosfera e gli elementi della sottostante superficie (suolo, acque, vegetazione) ne sono soggetti importanti.”

La definizione lascia nel vago i fenomeni micro-meteorologici, ma qui possiamo far ricorso all’esperienza di tutti i giorni:

  • E’ di sicuro appannaggio della micro-meteorologia la diffusione turbolenta, con applicazioni alla dispersione di inquinanti in atmosfera, ma anche di pollini, spore, animaletti.
  • Certamente, anche, appartengono alla micro-meteorologia i flussi turbolenti verticali (quindi tra superficie e atmosfera) di calore, vapore acqueo, anidride carbonica, metano, ammoniaca, …

Una presentazione divulgativa dei principali fenomeni micro-meteorologici, scritta in Italiano, si può trovare nel primo capitolo  del libro (Sozzi et al, 2002). Di questo bel testo, ormai da tempo fuori catalogo, è possibile trovare versioni PDF scaricabili, per esempio al sito di ARPA Lazio.

E quindi ci siamo…

Ciò che avviene sopra il tetto di Servizi Territorio dunque, appartiene a tutti gli effetti al campo della micro-meteorologia.

Di più: siamo certamente in città, e dunque in un ambiente micro-meteorologico di tipo urbano – in particolare, dei sotto-tipi “LCZ 2/3/E” (Oke et al, 2017).

Cos’è che compone una stazione micro-meteorologica?

Le componenti chiave di una stazione micro-meteorologica sono:

  • Un anemometro ultrasonico tri-assiale (come l’USA-1 della stazione di Cinisello Balsamo).
  • Un sistema di acquisizione ed elaborazione in tempo reale (una unità eddy-covariance MeteoFlux Core V2, nel caso della stazione di Cinisello).

A queste possono aggiungersene altre, ad esempio radiometri ed altri sensori più o meno esoterici.

Un primo esempio di dati (senza commento, per ora)

Le stazioni micro-meteorologiche misurano grandezze “in più”, rispetto alle stazioni meteorologiche tradizionali.

Una di queste grandezze in più, di cui in questa prima puntata non diciamo nulla, è la velocità verticale del vento.

Distribuzione empirica della componente verticale del vento.

Nelle puntate che seguiranno riprenderemo questa ed altre grandezze. E, cercheremo di trarne vantaggio.

Fonti e riferimenti

  • (Glickman, 2000) T.S. Glickman (ed.), Glossary of Meteorology, American Meteorological Society, 2000
  • (Oke et al., 2017) T.R. Oke, G. Mills, A. Christen, J.A. Voogt, Urban Climates, Cambridge University Press, 2017
  • (Sozzi et al., 2002) R. Sozzi, T. Georgiadis, M. Valentini, Introduzione alla turbolenza atmosferica: Concetti, stime, misure, Pitagora Editrice, 2002

 

I modelli di Lotka-Volterra generalizzati: storia di un fallimento fastidioso, ma molto, molto fecondo

Antefatto

La micro-meteorologia è una disciplina di confine, e come tutte le discipline di confine ha molte radici.

In questo articolo ne vedremo una, un po’ sorprendente per chi abbia un passato fisico. Ma, importante, e molto interessante – anche per quanto riguarda l’aspetto generale dell’utilità dei modelli matematici.

Una radice, biologica.

O più precisamente, sprofondata nel terreno dell’ecologia delle comunità.

Siamo circa alla fine degli anni ’80 del secolo scorso.

E’ un periodo ribollente, grazie alla diffusione dei primi calcolatori personali. Il calcolo aiuta la fantasia, e permette di immaginare applicazioni e tendenze che, prima, non erano neanche lontanamente concepibili.

E’, anche, un’epoca di grande ottimismo, almeno nel campo della biologia.  I progressi recenti nelle tecniche di osservazione hanno permesso di iniziare il sequenziamento dei geni più piccoli, e sulla scorta di questo progresso va consolidandosi un paradigma nuovo: la vita, come fenomeno deterministico. Quantificabile. Riducibile ad una descrizione algoritmica – e il meccanismo di traduzione dell’informazione contenuta nel DNA in proteine richiama in effetto molto da vicino un processo di codifica svolto da una Macchina di Turing.

Se una cosa del genere accadeva alla scala piccola degli individui, non poteva qualcosa di simile essere vero più in grande? Al livello di intere popolazioni, o di comunità ecologiche?

Col senno di poi, possiamo dire che tutto quell’ottimismo era in larga misura mal posto. Ma alcuni dei suoi strascichi hanno portato a scoperte molto interessanti. Tra queste, l’argomento di cui parliamo, i Modelli di Lotka-Volterra Generalizzati.

Predatori e prede: le equazioni di Lotka-Volterra

I Modelli di Lotka-Volterra Generalizzati sono una forma estesa del modello originale di Lotka-Volterra, introdotto indipendentemente nel corso degli anni Venti dello scorso secolo da Albert J. Lotka e Vito Volterra.

Il modello di Lotka-Volterra descrive l’interazione tra due popolazioni di specie diverse, in cui una si nutre a spese dell’altra.

Come in tutti i modelli di dinamica delle popolazioni in uso a quel tempo, anche il modello di Lotka-Volterra si esprime con un sistema di equazioni differenziali ordinarie del primo ordine:

In questa formula, i simboli Pp rappresentano rispettivamente le abbondanze della preda e del predatore; r è il tasso di accrescimento della popolazione preda, fosse lasciata in pace dal predatore; -s, un numero negativo, rappresenta il tasso di decrescita del predatore, decidesse di non nutrirsi a spese della preda. E infine, b è proporzionale all’entità di popolazione che la preda subisce ad opera della predazione, mentre kb è proporzionale al guadagno che la popolazione predatrice ricava per effetto della predazione.

Se facciamo finta che la preda e il predatore abbiano circa le stesse dimensioni e peso, è più che ragionevole assumere che il valore di k, una specie di rendimento di conversione della ciccia-di-preda in ciccia-di-predatore, sia positivo, ma minore di 1: stime condotte più o meno a spanne indicano per k un valore sensato dell’ordine di 0.1 – cioè, di tutta la biomassa che il predatore preleva dalla sua preda, solo un decimo si trasforma in biomassa del predatore.

Il bello del modello di Lotka-Volterra è che, in situazioni molto semplici, sembra funzionare. E spiega anche molto bene, in modo diretto, per un puro effetto matematico, l’andamento oscillante che si osserva in numerose coppie obbligate preda-predatore, specie in ambienti chiusi come nel caso di volpi e conigli nella taiga artica.

Il modello di Lotka-Volterra ha trovato anche altre applicazioni, per esempio nel campo dell’economia, con l’avvertenza di sostituire l’abbondanza di individui, o la biomassa, con il denaro, le prede con consumatori, e i predatori con produttori.

Il modello di Lotka-Volterra generalizzato

Ora: per quale motivo soltanto due popolazioni, e non di più? Perché, allora, non un’intera comunità trofica?

E perché, già che ci siamo, non ammettiamo la possibilità che la crescita dei produttori primari sia limitata da qualche fattore, come ad esempio lo spazio, o la quantità di azoto fissata nel terreno?

Se immaginiamo di abbracciare queste nuove ipotesi, arriviamo ad un’estensione del modello originale di Lotka-Volterra tutto sommato naturale e semplice da scrivere:

Questa formula è il “Modello di Lotka-Volterra Generalizzato” (GLV per gli amici). Le differenze rispetto alla versione originale sono:

  • Le specie possono essere più di due.
  • I tassi di crescita o decrescita, uno per ciascuna delle popolazioni, si possono immaginare arrangiati in un vettore, r.
  • Le interazioni tra le specie sono descritte da una “matrice di comunità”, C, quadrata, di ordine pari al numero di specie.
  • Gli elementi della matrice di comunità che si trovano sulla diagonale principale, e che descrivono l’effetto di feedback che una popolazione ha su sé stessa – cioè l’effetto dei fattori limitanti, possono essere diversi da zero (nel modello originale di Lotka-Volterra sono sempre nulli).
  • Gli elementi della matrice di comunità che si trovano su lati opposti della diagonale principale, e che nel caso del modello originale esprimevano il trasferimento (e la dissipazione) di biomassa dalla preda al predatore, adesso possono assumere qualsiasi valore. In particolare possono essere entrambi positivi (entrambe le specie ci guadagnano dall’interazione, come in una simbiosi mutualistica), o negativi (le due specie si combattono attivamente a vicenda, come accadeva tra i nostri antenati ed i grandi predatori della savana), oppure di segno opposto (una specie ci perde e l’altra ci guadagna, come nelle interazioni preda-predatore o ospite-parassita), oppure ancora, uno o due dei valori sono nulli (una delle specie non ha alcun effetto sull’altra).

Promettente, vero?

E, anche, molto facile da programmare con un calcolatore come quelli che si rendevano sempre più disponibili nella seconda metà degli anni Ottanta, con l’avvento del personal computer, e la sua diffusione nei dipartimenti universitari e nelle case private.

(Per inciso: i miei primi modelli di Lotka-Volterra generalizzati li ho codificati in una calcolatrice programmabile TI-58C, durante il mio periodo di internato volontario nel gruppo di ricerca del professor Guido Pacchetti, al Dipartimento di Biologia dell’Università Statale di Milano – era, più o meno, il 1983; la velocità di esecuzione non era proprio grandissima – ogni passo temporale richiedeva una quindicina di secondi – ma era comunque enormemente maggiore di quella che avrei potuto raggiungere facendo i calcoli a mano, con l’aiuto di una normale calcolatrice da tavolo, che poi ai miei tempi era ancora il metodo “normale”.)

Il problema della stabilità

Errori numerici a parte (per risolvere problemi differenziali come quelli dati dal modello di Lotka-Volterra generalizzato il “brutale” metodo di Eulero non è il massimo, e metodi tipo Rune-Kutta sono decisamente preferibili – ma, con le risorse allora disponibili, non si sottilizzava troppo), i modelli di Lotka-Volterra generalizzati (GLV per gli amici) mostravano un piccolo, trascurabile difetto – diciamo così. Se provavi ad assegnare la matrice di comunità C nel modo più realistico possibile, e attribuivi alle abbondanze delle specie un valore iniziale inversamente (più o meno) proporzionale al livello nella catena alimentare, molto spesso dopo un numero anche piccolo di passi temporali una o più delle specie assumevano un’abbondanza negativa: si estinguevano.

Magari, nel frattempo, qualche altra specie divergeva e cominciava a crescere esponenzialmente (cosa che, a mia memoria, accadeva piuttosto di rado, e sempre per qualche errore nell’impostazione della matrice di comunità o del vettore dei tassi di variazione intrinseca).

Si vede che i “nostri” problemi non erano solo conseguenza di qualche errore sistematico, o di un effetto dovuto all’aria milanese: piuttosto presto comparvero articoli nei quali altri gruppi, piuttosto allarmati, dichiaravano problemi di stabilità.

Strano, data la ricchezza espressiva del modello GLV, che pareva riuscisse a racchiudere tutte le possibili interazioni tra specie e specie, e che a priori avrebbe dovuto poter descrivere facilmente la dinamica di ogni ecosistema.

Ci furono, ricordo, vari tentativi di miglioramento. Per esempio, passando dalla rappresentazione delle popolazioni per abbondanza a quella per biomassa. O, riscrivendo le equazioni applicando astuti metodi di scaling.

Ma, niente. Le comunità simulate dai modelli GLV continuavano a perdere pezzi (specie) per strada, sino a ridursi al nulla. Le poche che davano qualche segno di stabilità erano artificiali, semplicissime, con modalità di trasferimento di biomassa peculiari.

E non basta: alcune combinazioni di matrice di comunità e tassi di variazione intrinseci davano luogo a soluzioni caotiche – molto belle e stimolanti sul piano matematico, ma non proprio esattamente vicinissime  a quello che si può vedere in Natura.

Era piuttosto ovvio che qualcosa sfuggiva. Qualcosa di macroscopico.

Muore un’idea vecchia. Nasce un’idea nuova.

Proprio nel momento di massima perplessità (era divenuto sin troppo evidente che i modelli GLV non spiegavano la complessità degli ecosistemi che si trovano su questo pianeta), lasciai il gruppo del professor Pacchetti – dovevo concentrarmi sui miei studi, cercando con fatica di prendere la mia laurea in matematica, indirizzo applicativo, e già che c’ero cercarmi un lavoro.

Ma non persi del tutto i contatti, e potei raccogliere la storia di “com’è andata dopo” direttamente dalla voce di Guido, e, indirettamente, dalla storia delle tesi di alcuni suoi laureandi.

In effetti, il “baco macroscopico” c’era. Lì, in piena luce, ma nessuno aveva avuto modo di pensarci, tutti i gruppi essendo presi dalle relazioni tra specie.

Il problemone, in effetti, era l’interazione delle specie con l’ambiente fisico.

Al tempo dello sviluppo dei modelli GLV, con un atto molto antropomorfico, si buttava tutto il cuore (anche la ricerca scientifica lo ha, non di rado piuttosto conservatore, o quantomeno condizionato da sorprendenti elementi di soggettività) dalla parte dei predatori. Perché, primo, noi umani siamo il Predatore di Vertice; e, secondo, tutti i predatori non umani ci ispirano un senso di trasporto, di empatia, che difficilmente riserveremmo a dei brucatori, o a delle piante. Non per caso, i nostri amici a quattro zampe sono cani e gatti – predatori, con una testa da predatore, e quindi con tanto in comune con noi – mentre le mucche le alleviamo (per mangiarle), e le piante le coltiviamo (idem).

Noi predatori, vuoto per pieno, ricaviamo tutto il nostro sostentamento dalle nostre prede naturali. E una comunità trofica “a misura di testa di predatore” sarebbe una rete di relazioni tra specie in cui ognuna ricava il suo sostentamento dalle altre. Insomma: sarebbe descritta da un modello GLV, e “chiuso il cinema”.

Ma non ci sono solo i predatori.

Ci sono anche le prede.

E, tra queste, quelle alla base della base di tutta la catena alimentare: le piante verdi, i cianobatteri.

Le creature autotrofe.

Che da qualche parte devono pure ricavarlo, il loro sostentamento.

Ed in effetti lo ricavano, scambiando massa ed energia con il “mondo minerale”.

In effetti, gli scambi tra mondo minerale e biosfera sono ingentissimi, almeno un ordine di grandezza maggiori di tutti gli scambi di tutte le specie tra loro.

Questi scambi, di enorme dimensione, avvengono con tranquilla continuità, con molta calma, e noi frenetici animali non siamo nella posizione migliore per percepirli. Avvengono troppo “lentamente”.

Ma è grazie ad essi, che le piante verdi e i cianobatteri regolano attivamente la percentuale di ossigeno nell’atmosfera, a loro uso, per crearsi un ambiente ottimale nel quale la fotosintesi clorofilliana può procedere senza intoppi, e non c’è il pericolo di una sua inversione in “fotorespirazione”, cosa che accadrebbe ad una concentrazione di ossigeno troppo alta. Oppure, nel quale l’acqua fluisce dalle piante verso l’atmosfera (rendendo conto di circa il 70% dell’evaporazione sulla terraferma) ad un rateo compatibile con le necessità delle cellule vegetali.

Insomma: i modelli GLV davano per scontata una sorta di posizione dominante da parte di noi umani e delle creature più simili a noi, nell’economia generale della vita. La realtà, però, è molto diversa: la fama di vita davvero dominante è quella “vegetale” in senso lato, e noi animali siamo un piccolo disturbo (funzionale alla sopravvivenza delle piante, e soprattutto alla loro decomposizione). Piante, cianobatteri, archeobatteri, a loro volta sono impegnati in colossali scambi nei quali giocano un ruolo determinante numerosi fenomeni fisici.

Ad esempio, la turbolenza negli strati più bassi dell’atmosfera.

Questa constatazione, originata dal fallimento dei modelli GLV, ha avuto conseguenze di grandissima importanza, innescando una sorta di rivoluzione copernicana.

Flussi di (materia ed) energia

Quindi: la descrizione in termini di scambio di biomassa tra specie viventi non funzionava.

Ma una descrizione più completa, in termini di flussi di energia (e di un po’ di biomassa, che in fondo è una sorta di energia immagazzinata in forma chimica, e di sostanze inorganiche varie), poteva invece dare risposte più sensate.

E così, in effetti, è accaduto: gli ecosistemi hanno trovato una descrizione più efficace, in termini di flussi di energia e di materia.

La matematica che sta dietro a questa descrizione, inutile dire, è più complicata di quella dell’algebra lineare dei modelli GLV: adesso, si usano sistemi di equazioni differenziali alle derivate parziali.

E’, anche, nata una nuova disciplina, l’ecologia degli ecosistemi, che descrive l’interazione tra biosfera e sistema-Terra visto nel suo complesso, con variabili indipendenti che non sono più (soltanto) abbondanze di specie o biomasse, ma, grandi flussi di nutrienti ed energia a scala globale, o di grandi regioni geografiche.

La necessità di studiare i flussi di massa e di energia che occorrono tra grandi sotto-sistemi, nel caso tra superficie della Terra e biosfera da un lato, ed atmosfera dall’altro, ha contribuito alla fioritura di interesse per tecniche di misura come la eddy covariance, che oggi trova applicazioni non solo nel campo della fisica dell’atmosfera, ma anche nelle scienze agricole e forestali, con un’attenzione particolare ai flussi di calore, vapore acqueo, anidride carbonica e, ultimamente, anche metano ed ammoniaca.

Tutti gli studi in corso contribuiscono a darci un quadro del sistema-vita molto più articolato di quanto accadeva al tempo dei modelli GLV. Un quadro, nel quale le connessioni interessano non soltanto l’interno della biosfera, ma elementi spesso sorprendenti del sistema-Terra, del Sole, dell’atmosfera e dell’idrosfera. Su scale temporali che, ormai, coprono un intervallo che va dalle centinaia di milioni di anni al decimo di secondo.

Molto c’è ancora da scoprire, di misterioso ed affascinante. In questi ultimi anni, ad esempio, stanno moltiplicandosi studi nei quali si cerca di indagare il modo, apparentemente molto raffinato, in cui le molecole della vita riescono a sfruttare fenomeni di natura quantistica nonostante le loro grandi dimensioni, ottimizzando reazioni (come la fotosintesi) che, in un mondo puramente classico, non sarebbero così singolarmente efficienti come sono invece in Natura.

Insomma, stiamo comprendendo sempre più di non sapere.

Una cosa, comunque, risulta chiara: i modelli GLV erano senza dubbio “sbagliati”, ma il loro errore è stato, come minimo, felice. Il superarlo, ha aperto strade che non potevamo, prima, immaginare.

L’uso delle riprese time-lapse delle nubi come mezzo per visualizzare la circolazione locale (e non) dell’aria

 Una ricerca superficiale su Internet ci può fare facilmente scoprire centinaia di webcam disseminate sul territorio.

Molte di queste non si limitano ad acquisire immagini statiche di tanto in tanto, ma sono anche in grado di registrare filmati, spesso nella modalità “time lapse“.

Per chi ancora non ne abbia veduto uno, qui sotto ho riportato un esempio, un collage di tre inquadrature riprese da Moggio Valsassina guardando verso le Grigne ed i Piani di Artavaggio. In pratica, un filmato time lapse è costituito dall’assemblaggio di fotografie scattate ad intervalli regolari nel tempo.

La relazione tra il rateo di proiezione delle fotografie nel filmato e l’intervallo di tempo reale tra foto consecutive determina un effetto di accelerazione.

Il filmato qui sopra, per esempio, è stato ottenuto scattando le singole fotografie ogni 2 secondi, con un rateo di proiezione di 30 immagini al secondo. Questo vuol dire che un secondo di tempo di proiezione (30 fotogrammi) equivale a 30*2, cioè 60 secondi di tempo reale, per un’accelerazione apparente del tempo pari a 60 volte.

L’accelerazione apparente del tempo permette di seguire l’evoluzione di fenomeni troppo lenti perché noi possiamo percepirli ad occhio nudo, come ad esempio l’aprirsi e chiudersi dei fiori, il crescere dei vegetali, o lo svilupparsi di una grande struttura in costruzione, come un edificio od una nave.

Oppure, anche, il muoversi ed il continuo trasformarsi delle nubi nel cielo.

In moltissimi casi, il movimento apparente delle nubi risulta dalla complessa interazione tra il flusso dell’aria, con il suo complicato apparire a quote diverse ed in presenza di ostacoli naturali come le montagne, e i delicati equilibri dinamici dovuti all’alternanza di condensazione ed evaporazione – come accade per esempio con le nubi orografiche che sembrano stazionare sulle cime dei rilievi nelle giornate calde di bel tempo, e che sono dovuti al continuo formarsi e dissolversi di nubi da condensazione, mentre l’aria investe la montagna, si innalza raffreddandosi, per poi ridiscendere riscaldandosi, con effetto sulla condensazione.

Tutti questi effetti, ad ogni modo, sono chiaramente riconoscibili e distinguibili in un filmato time lapse. Già nell’esempio che ho riportato possiamo distinguere facilmente regimi di circolazione diversi alle varie quote. Figuriamoci cosa potrebbe coglierne un occhio esperto, e quale appoggio potrebbe dare, con immagini che come spesso accade valgono centinaia di parole, ad una conferenza, o una lezione.

Certo, perché le immagini possano avere un reale uso didattico è bene soddisfino ad alcuni requisiti di base – che non è detto i filmati raccolti automaticamente dalle webcam abbiano in pieno.

Per esempio, esposizione e contrasto non dovrebbero mostrare eccessivi “salti” durante la proiezione. Poi, l’immagine dovrebbe fluire con regolarità, e non a scatti. In presenza di aree con luminosità molto diverse, le immagini non dovrebbero oscurarne alcune e bruciarne altre.

Insomma: per ottenere risultati utili serve un po’ di lavoro fotografico.

Ma la tecnologia di oggi viene in soccorso, con programmi di assemblaggio ed editing delle sequenze time lapse, o addirittura con videocamere dedicate (a basso costo, oltretutto – tra le webcam di fascia media e le macchine fotografiche compatte). Non desidero, in questa sede, fare pubblicità agli uni od agli altri, ma una ricerca anche breve su Internet vi mostrerà decine di opportunità.

Certo, che vedere immagini di nubi in movimento fa senz’altro pensare e, allo stesso tempo, apprezzare (in modo molto rilassante) la meraviglia del mondo.