Introduction
In France, main geological/geoscientific data are managed by the French Geological Survey and several platforms provide a public and free access to many datasets such as:
- geological maps,
- boreholes, piezometers and wells locations, geological reports and cross sections (InfoTerre). In France, this database is known as the BSS (Base de donnée du Sous-Sol) which includes more than 700 000 boreholes descriptions,
- piezometric time series and physico-chemical analysis from the national monitoring network (ADES).
These data dedicated to the (hydro-)geological knowledge are used by engineers and scientists to prevent natural risks and disasters (e.g. landslides or settlements), protect groundwater ressources and secure water supply for example. Particularly, these data can be used in many fields such as hydrogeology, geothermal energy and geotechnical applications.
These datasets can be explored and download manually by visiting the main platform of the french geological survey InfoTerre, but when it comes to export a large amount of data, at a regional scale, this task becomes tedious and time consuming. To face this problem and improve the dissemination of geoscientific data to the public and to scientists, the French geological survey recently developed many geo-webservices in accordance with European open-data requirements. These geo-webservices can help to automate geological data manipulation, analysis and mapping.
Consequently, the aim of this article is to present some of these services. Particularly, we are going to:
- download and analyze boreholes information at the scale of the department of our choice: water supply wells and geothermal installations will be identified.
- map all this hydrogeological information on an interactive folium map and add the geological tiles of France.
- play with two Hubeau APIs to describe piezometric and physico-chemical characteristics of some groundwater wells: we will learn how to access and analyze hydraulic heads and groundwater quality fluctuations on the groundwater wells we want.
Run me first
Install/import usefull libraries
Before starting we need to install/import some libraries:
Administrative units of France
To download information at the scale of a department, we need some a geographical dataset with administrative units of France. We can use the openstreetmap shapefile of departments of France. We download this shapefile and put it in a local directory as follow:
code_insee | nom | nuts3 | wikipedia | surf_km2 | geometry | |
0 | 974 | La Réunion | FR940 | fr:La Réunion | 2505 | MULTIPOLYGON (((55.21643 -21.03904, 55.21652 -… |
1 | 11 | Aude | FR811 | fr:Aude (département) | 6343 | POLYGON ((1.68872 43.27368, 1.69001 43.27423, … |
2 | 43 | Haute-Loire | FR723 | fr:Haute-Loire | 5003 | POLYGON ((3.08206 45.28988, 3.08209 45.29031, … |
3 | 13 | Bouches-du-Rhône | FR823 | fr:Bouches-du-Rhône | 5247 | MULTIPOLYGON (((4.23014 43.46047, 4.23025 43.4… |
4 | 47 | Lot-et-Garonne | FR614 | fr:Lot-et-Garonne | 5385 | POLYGON ((-0.14058 44.22648, -0.12931 44.23218… |
All department numbers are indicated in the column code_insee. We should be careful about this number because the Rhône department is identified as 69D to make the difference with the metropolitan area of Lyon which is identified as 69M. Here, we do not want to distinguish them so we create a new column associating the same number 69.
Selection of some departments of interest
To see the result, we can select some departments of interest and preview their shape. In the following, we will work with this selection of departments numbers to get some (hydro-)geological information. I decided to work with the department of the Ain (01), the Rhône (69) and of the Isère (38).
We also download a dataset of main cities of France (more than 100’000 inhabitants) and we keep the cities inside our depatment of interest following a similar procedure.
city | lat | lng | country | iso2 | admin_name | capital | population | population_proper | geometry | index_right |
8 | Lyon | 45.76 | 4.84 | France | FR | Auvergne-Rhône-Alpes | admin | 516092 | 516092 | POINT (4.84000 45.76000) |
5 | Grenoble | 45.1715 | 5.7224 | France | FR | Auvergne-Rhône-Alpes | minor | 687985 | 158454 | POINT (5.72240 45.17150) |
165 | Bourg-en-Bresse | 46.2056 | 5.2289 | France | FR | Auvergne-Rhône-Alpes | minor | 41527 | 41527 | POINT (5.22890 46.20560) |
256 | Vienne | 45.5242 | 4.8781 | France | FR | Auvergne-Rhône-Alpes | minor | 29306 | 29306 | POINT (4.87810 45.52420) |
Extraction and analysis of boreholes information (BSS data)
Get our data into a clean (geo-)dataframe
Now we can make a loop and concatenate the data of our departments of interest:
ID_BSS | indice | designation | x_ref06 | y_ref06 | lex_projection_ref06 | … | code_bss | |
0 | BSS001PYDG | 06246X0004 | CPT | 810400 | 6575119 | Lambert-93 | … | 06246X0004/CPT |
1 | BSS001PYDH | 06246X0005 | CPT | 810316 | 6574955 | Lambert-93 | … | 06246X0005/CPT |
2 | BSS001PYDJ | 06246X0006 | CPT | 810209 | 6575025 | Lambert-93 | … | 06246X0006/CPT |
3 | BSS001PYER | 06248X0003 | S | 827507 | 6575474 | Lambert-93 | … | 06248X0003/S |
4 | BSS001PYES | 06248X0004 | HY | 828134 | 6576416 | Lambert-93 | … | 06248X0004/HY |
As you can see, there are 71 columns. In the following, we reduce our focus to few relevent attributs (but please add more if you want to):
- ID_BSS, indice and designation which are national codes of boreholes. These codes are requested when making queries to different webservices (we’ll see that later),
- point_eau is a bolean indicating if there is some groundwater or not,
- x_ref06 and y_ref06 are indicating the location in a geographical projection given by lex_projection_ref06,
- z_sol indicates the ground elevation mesured and z_bdalti the ground elevation given by the IGN at the same location (can be useful to get both to identify some errors),
- prof_investigation indicates the total depth of the borehole,
- prof_eau_sol indicates the groundwater table depth if any exists,
- lex_usage refers to the usage of the borehole/water well,
- procede_geothermique, categorie_geothermie, and usage_geothermie are indicating if the point is actually a geothermal installation.
Your content goes here. Edit or remove this text inline or in the module Content settings. You can also style every aspect of this content in the module Design settings and even apply custom CSS to this text in the module Advanced settings.
Because the projection is not the same for all departments of France (because of overseas departments), the epsg number of different projections must be known. They are listed in a epsg dataframe below and finally, a procedure is made to convert all location in longitudes/latitudes (epsg:4326). The following step is required in case you want to work on overseas departments. In case you won’t, the following procedure still works to get longitudes/latitudes from Lambert 93 coordinates.

We can of course save our geodataframe as a shapefile, in a RESULTS folder as follow:
Identifying groundwater wells and geothermal installations
Of course we can decide to distinguish the kind of borehole we want to keep. For example, let’s select (1) boreholes where a groundwater table depth is reported and (2) geothermal installations.
In the first case, we keep points where the attribute point_eau is OUI (meaning that there is some groundwater right here) AND where prof_eau_sol is not a NaN value (meaning that we should know the groundwater table depth). In the second case, we keep points where categorie_geothermie or procede_geothermique or usage_geothermie are not a NaN value.

Describing groundwater table elevation
Now that we know where groundwater bodies are located, we might want to describe groundwater levels and maybe describe groundwater flow directions. To do so, we refine a location of interest by choosing the coordinate of a point, we also define a relevant buffer zone around this point to reduce our dataset to this area of interest. In our example, we select the city of Lyon, and a buffer zone of 10km.
Please note that groundwater table elevation fluctuates over time. In the following, the calculated value represents the water table elevation observed when the borehole was drilled. Also, it should be noted that the calculated values are not synchronous which can alter the local representativity regarding the interpretation groundwater flow directions and flow systems.
Please note that we could instead import the polygon shapefile of our choice and select boreholes inside our polygon.
We describe data inside our geodataframe:
z_sol | z_bdalti | prof_investigation | prof_eau_sol | |
count | 2080 | 1860 | 2008 | 2092 |
mean | 174.73377 | 181.92688 | 19.682102 | 8.022075 |
std | 86.678638 | 30.315466 | 13.985364 | 6.130222 |
min | -999 | 155 | 2.3 | 0.2 |
25% | 166.7975 | 167 | 12 | 4.2575 |
50% | 169.92 | 170 | 17.5 | 6.1 |
75% | 182 | 184 | 22.5 | 10 |
max | 370 | 379 | 250.25 | 85 |
The result shows some errors regarding the ground elevation (-999 value). Consequently, let’s remove data outside first and last percentiles:
z_sol | z_bdalti | prof_investigation | prof_eau_sol | |
count | 2036 | 1826 | 1953 | 2036 |
mean | 178.90788 | 180.04272 | 19.775622 | 8.048355 |
std | 23.567719 | 24.353793 | 14.11073 | 6.157885 |
min | 157.044 | 155 | 2.3 | 0.2 |
25% | 166.9 | 167 | 12 | 4.3 |
50% | 169.945 | 170 | 17.5 | 6.1 |
75% | 182 | 183 | 23 | 10 |
max | 329 | 330 | 250.25 | 85 |
Now, we define the groundwater table elevation as the ground elevation minus the groundwater table depth. Because we have two different possibility to choose the ground elevation from the IGN z_bdalti. In addition, we can consider that in our area of interest we focus on the shallow aquifer where groundwater is flowing in the modern alluvial deposit. Consequently, let’s focus on an investigation depth ranging from 10m to 30m:
Then, groundwater table elevation values can be displayed as follow:

Interactive mapping with folium
The information above can easily be projected on an interactive map using the folium library. Additionally, it is interesting to do so because it can help us to understand the relation with geological units.
In fact, since a few years, geological tiles are available across France and can be projected and visualized using folium. All WMF/WFS services provided by the BRGM are listed in a table at this page. It is indicated that the geological map can be visualized using the following URL http://geoservices.brgm.fr/geologie, and the following layer GEOLOGIE. We directly use these inputs to define our WMS tile layer (see below in the script).
Then we can build our map using folium:
Of couse we can save this map in html and open it later with any navigator:
Access to hydraulic head monitoring data
Since a couple of years, several APIs dedicated to water management have been developed on a national platform called Hubeau. Some of these APIs are dedicated to hydrobiology and surface water monitoring, and others to quantitative and qualitative description of groundwater.
The principle is quite easy: you provide the API of your choice with (1) one or several borehole identifiers, (2) a period of interest, (3) some parameters of interest (piezometric level, temperature, etc.) and you get back a dictionary with all information.
First the procedure to describe the hydraulic head fluctuation at a given borehole will be described. Secondly, we will focus on how to find all recording stations around a area of interest, for a given period of time.
Getting hydraulic head fluctuations
Let’s see how the groundwater piezometric API works. We first need to define one borehole identifier code_bss, and a period of interest:
To query the API, we need to build an URL which depends on previous parameters. I give the function below providing you the appropriate URL:
We can now call this function as follow:
It returns the following URL you can visit. We can now make our query using the resquests.get function:
The code ‘200’ means that our query is successful (the code is valid and there are some data available). We now organize the result (using the json method) and print the content of the webpage using the following procedure:
The resulting content is a dictionary organized according to the date of measurements. We now need to define a function to extract the information we want in this dictionary. The first key of this dictionary gives us a count of all measurements over our period of interest. All piezometric values are stored in a data sub-dictionary. Let’s see on the first element of the data dictionary:
{‘code_bss’: ‘06987A0186/S’,
‘code_continuite’: ‘2’,
‘code_nature_mesure’: None,
‘code_producteur’: ‘196’,
‘date_mesure’: ‘2015-01-01’,
‘mode_obtention’: ‘Valeur mesurée’,
‘niveau_nappe_eau’: 163.09,
‘nom_continuite’: ‘Point lié au point précédent’,
‘nom_nature_mesure’: None,
‘nom_producteur’: ‘Service Géologique Régional Rhône-Alpes (196)’,
‘profondeur_nappe’: 5.41,
‘qualification’: ‘Correcte’,
‘statut’: ‘Donnée contrôlée niveau 2’,
‘timestamp_mesure’: 1420117200000,
‘urn_bss’: ‘http://services.ades.eaufrance.fr/pointeau/06987A0186/S‘}
We have here a new dictionary with the ID of our borehole, the date of the measurement (date_mesure), the associated timestamp (timestamp_mesure), the hydraulic head measured (niveau_nappe_eau) and other information.
We can make a function to convert this dictionary into a dataframe keeping only relevant features as follow:
date | hydraulic_head |
01/01/2015 | 163.09 |
02/01/2015 | 163.08 |
03/01/2015 | 163.07 |
04/01/2015 | 163.06 |
05/01/2015 | 163.19 |
Finding piezometric monitoring stations around a location of interest
A procedure exists to find all stations with hydraulic head records around a location of interest. What you need is to define the longitude and latitude coordinates of a rectangle area. Let’s take the smallest rectangle including our previous circular buffer zone:
To make our query, we only need the lower left and top right corners coordinates:
We then need to make a new function to produce the URL and another one to organize the result in a pandas dataframe:
We run this function to sew how many piezometric stations are working/recording at a given date: the 1st of january 2020.
code_bss | date_debut_mesure | date_fin_mesure | nom_commune | x | y | bss_id |
06988C0281/F | 10/04/2009 | 17/08/2020 | Chassieu | 4.969186 | 45.742375 | BSS001TPVW |
07223X0130/P | 18/05/2005 | 03/03/2020 | Saint-Priest | 4.908662 | 45.701409 | BSS001URVP |
06987X0272/P | 14/09/2005 | 09/03/2020 | Vaulx-en-Velin | 4.937039 | 45.807487 | BSS001TPEW |
06987X0274/P | 14/09/2005 | 09/03/2020 | Vaulx-en-Velin | 4.910798 | 45.796308 | BSS001TPEY |
06988X0217/P | 23/08/2005 | 09/03/2020 | Décines-Charpieu | 4.952973 | 45.776655 | BSS001TQGK |
06987A0186/S | 16/04/1968 | 07/08/2020 | Villeurbanne | 4.864856 | 45.779492 | BSS001TMCR |
06987J0105/PZ | 16/05/2005 | 03/03/2020 | Villeurbanne | 4.899866 | 45.782911 | BSS001TNEM |
We can of course display the location of these active stations on an interactive folium map:
Access to physico-chemical parameters
Similarly, we can make some queries to observe the fluctuation of some physico-chemical parameters such as groundwater temperature, nitrates concentrations, pesticides… many other parameters are available.
Exploration of accessible physico-chemical parameters
Since 1992, the French producers of public water data have engaged a consistency process for their data in the framework of the French water information system. Consequently, all groundwater quality parameters are associated to a number in a referential SANDRE.
This referential can be download below. We can see that 1845 parameters are available.
Nom_famille_sise | nom_parametre_sise | code_parametre_sandre | |
1 | COMPOSES ORGANOHALOGENES VOLATILS | Tétrachloroéthane-1,1,1,2 | 1270 |
2 | COMPOSES ORGANOHALOGENES VOLATILS | 1,1,1,2 Tétrachloropropane | 2704 |
3 | COMPOSES ORGANOHALOGENES VOLATILS | 1,1,1,3 Tétrachloropropane | 2705 |
4 | COMPOSES ORGANOHALOGENES VOLATILS | Trichloroéthane-1,1,1 | 1284 |
5 | COMPOSES ORGANOHALOGENES VOLATILS | Tétrachloroéthane-1,1,2,2 | 1271 |
… | … | … | … |
1841 | PHYTOPLANCTONS | % de colo de ulothricophyc. subst | 6431 |
1842 | PHYTOPLANCTONS | Colonies de zygophyc. substrat | 6432 |
1843 | PHYTOPLANCTONS | % de colo de zygophyc. substr | 6432 |
1844 | PARAMETRES MICROBIOLOGIQUES | Entérocoques / 100 mL (qPCR) | 6455 |
1845 | PHYTOPLANCTONS | Radiocystis sp (cellules) substrat | 7299 |
We could explore all accessible parameters, but in this tutorial, we focus on two of them:
- temperature with the code ‘1301’,
- nitrates concentration with the code ‘1340’.
It can be checked by printing the appropriate lines of the dataframe:
Nom_famille_sise | nom_parametre_sise | code_parametre_sandre | |
1255 | PARAMETRES AZOTES ET PHOSPHORES | Nitrates (en NO3) | 1340 |
1645 | CONTEXTE ENVIRONNEMENTAL | Température de l’eau | 1301 |
The information about these parameters (units, etc.) are described in the following webpages (in french):
You can easly access and get the information about the parameter you want by using the following function depending on the parameter code:
Finding quality monitoring stations around a location of interest
First, we find some quality monitoring stations around our location of interest by using the function defined earlier:
bss_id | code_bss | date_debut_mesure | date_fin_mesure | longitude | latitude | altitude | nom_commune | … | |
0 | BSS001TLEA | 06986T0059/S7 | 05/02/2004 | 05/02/2004 | 4.840702 | 45.731292 | 168 | Lyon 7e Arrondissement | … |
1 | BSS001TPEA | 06987X0252/PZ4 | 14/10/2002 | 02/06/2004 | 4.875022 | 45.762547 | 168 | Villeurbanne | … |
2 | BSS001TPFA | 06987X0276/PZ2 | 31/10/2003 | 04/10/2004 | 4.869878 | 45.750762 | 170 | Lyon | … |
3 | BSS001URFA | 07222X0490/PN4 | 27/11/2000 | 12/12/2007 | 4.845675 | 45.692831 | 164 | Saint-Fons | … |
4 | BSS001URGA | 07222X0514/PZ1 | None | None | 4.845893 | 45.706085 | 164 | Saint-Fons | … |
We see that there are 471 stations. However, the dates are not necessarily matching our query (it is probably an issue that will be solved in a future version of the API). To refine our selection on our period of interest we convert all date features in datetimes format and compare with our datetimes of interest:
We have finally 7 active stations and we can add them (red color) to our interactive map:
Getting groundwater temperature and nitrates concentration fluctuations
We first need to define one borehole identifier code_bss, a period of interest and the code of parameters we want to explore (temperature and nitrates in our example):
To query the API, we need to build an URL which depends on previous parameters. I give the function below providing you the appropriate URL depending on the period and the parameters we want:
We can now call this function as follow:
Now, we get and organize the content of this webpage into a pandas dataframe. Please note that a very large number of features are available. You can explore these features by visiting the API website or by removing the last command in the following cell.
Here we only keep following features:
- the date of the measurement,
- the code of the parameter measured,
- the name of the parameter measured (in french),
- the resulting value,
- the quantification limit,
- the detection limit,
- the analytical uncertainty
and we finally print our dataset:
date | code_param | nom_param | value | limite_quantification | limite_detection | incertitude_analytique |
11/05/2017 | 1340 | Nitrates | 19 | 0.1 | None | None |
07/12/2017 | 1340 | Nitrates | 19 | 0.1 | None | None |
07/12/2017 | 1301 | Température de l’Eau | 16.7 | NaN | None | None |
05/07/2018 | 1340 | Nitrates | 17.5 | 0.1 | None | None |
05/07/2018 | 1301 | Température de l’Eau | 17 | NaN | None | None |
16/10/2018 | 1301 | Température de l’Eau | 17 | NaN | None | None |
16/10/2018 | 1340 | Nitrates | 17.9 | 0.1 | None | None |
21/05/2019 | 1340 | Nitrates | 20 | 0.25 | None | None |
21/05/2019 | 1301 | Température de l’Eau | 17 | NaN | None | None |
16/10/2019 | 1301 | Température de l’Eau | 17.2 | NaN | None | None |
16/10/2019 | 1340 | Nitrates | 18 | 0.25 | None | None |
Discussion
This notebook intends to introduce some APIs and geo-webservices dedicated to the (hydro-)geological knowledge of France. This notebook has been reviewed and successfully returns several parameter descriptions around an area of interest. Please note that all input data have been carefully checked, and be aware that an alteration of input data might raise errors and/or inconsistencies. Particularly, some site-specific choices have been made to illustrate different possibilities and must not be generalized to other locations. Of course, readers are invited to play with all input parameters, and focus on other locations, but following advises and comments are given:
- the area of interest should be carefully determined and must be of a reasonable size to limit the execution time of different procedures. Particularly, launching a data extraction of the all France is misadvised,
- the closeness and depth of boreholes must be taken carefully when comparing some measurements: un-correlated fluctuations can be observed when sensors are note located into the same groundwater body or geological unit, even if boreholes are close.
- there are many more available API functionalities not discussed in this notebook. Please visit Hubeau and PIC’EAU to find out more.
- to work at reginal scales on several stations at a time, there are more efficient methods to extract and analyse data than ones presented in this notebook.
Play with the script…
This post is available as a notebook and can be run on Goolge Colab.
Please note that all my notebooks are available on my GitHub repository PythonForGeosciences.