Map Volcanic Ash Cloud Advisories: Wolfram Language Code Gallery

Map Volcanic Ash Cloud Advisories

Make maps of the volcanic ash clouds reported by Volcanic Ash Advisory Centers.

code

MapVAAdvisory[msg_] := Module[{observationString, observation, observationLayers, polys, date}, observationString = StringCases[msg, "OBS VA CLD: " ~~ Shortest[x__] ~~ Alternatives["FCST VA CLD", EndOfString] :> x]; If[observationString === {}, Return[Null]]; observation = StringReplace[observationString[[1]], "\n" -> " "]; If[StringCases[observation, RegularExpression["[NESW][0-9]{4} [NESW][0-9]{5} "]] === {}, Return[Null]]; observationLayers = StringSplit[ observation, {"SFC/FL" ~~ DigitCharacter .. ~~ Whitespace, "FL" ~~ DigitCharacter .. ~~ "/" ~~ DigitCharacter .. ~~ Whitespace}]; polys = (GeoPosition@FromDMS[StringJoin[ StringTake[#, 3], "d", StringTake[#, {4, 5}], "m", StringTake[#, {6, 10}], "d", StringTake[#, {11, 12}], "m"]] & /@ StringSplit[#, " - "]) & /@ observationLayers; date = StringCases[msg, "DTG: " ~~ Shortest[x__] ~~ Alternatives[Whitespace, "\n"] :> x][[1]]; GeoGraphics[{ GeoStyling[Opacity[.5]], MapIndexed[With[{color = {Yellow, Red, Orange}[[#2[[1]]]]}, {EdgeForm[color], FaceForm[color], Polygon[#1]}] &, polys] }, PlotLabel -> "Volcanic Ash Cloud Over Europe\n" <> date] ];

how it works

The Volcanic Ash Advisory Centers (VAAC) issue advisories when civil aviation is endangered by volcanic eruptions. In their textual form, advisories are not especially enlightening, but you can parse the data in the advisory and make a map of the ash cloud using the pattern matching, string processing, and geomapping capabilities built into the Wolfram Language.

Here’s an example of an advisory from the Toulouse VAAC for an ash cloud resulting from an eruption of the Furnas volcano in the Azores:

msg = Import[ "http://www.meteo.fr/vaac/archives/FURNAS/e.20120427120300.FURNAS.\ 201207.201204271200.html"]

Extract the part of the advisory that specifies the observed cloud:

observationString = StringCases[msg, "OBS VA CLD: " ~~ Shortest[x__] ~~ Alternatives["FCST VA CLD", EndOfString] :> x]

The first element of that list is the observation; replace any new lines with spaces so that the string is uniformly formatted:

observation = StringReplace[observationString[[1]], "\n" -> " "]

The observation consists of cloud extents for multiple atmospheric layers. Split the observation into separate layers to obtain separate lists of latitude longitude pairs:

observationLayers = StringSplit[observation, {"SFC/FL" ~~ DigitCharacter .. ~~ Whitespace, "FL" ~~ DigitCharacter .. ~~ "/" ~~ DigitCharacter .. ~~ Whitespace}]

Process each layer string to obtain lists of {latitude, longitude} GeoPositions:

polys = (GeoPosition@FromDMS[StringJoin[ StringTake[#, 3], "d", StringTake[#, {4, 5}], "m", StringTake[#, {6, 10}], "d", StringTake[#, {11, 12}], "m"]] & /@ StringSplit[#, " - "]) & /@ observationLayers

Extract the date and time from the observation:

date = StringCases[msg, "DTG: " ~~ Shortest[x__] ~~ Alternatives[Whitespace, "\n"] :> x][[ 1]]

Use the extracted date and cloud extents to make a labeled map of the ash cloud:

GeoGraphics[{ GeoStyling[Opacity[.5]], MapIndexed[With[{color = {Yellow, Red, Orange}[[#2[[1]]]]}, {EdgeForm[color], FaceForm[color], Polygon[#1]}] &, polys] }, PlotLabel -> "Volcanic Ash Cloud Over Europe\n" <> date]