"NotAvailable" and "NotApplicable" data points returned by WeatherData[]

Question

I am trying to use ListVectorPlot in combination with WeatherData and CountryData to show wind velocity vectors over a given country. The problem I am having stems from "NotAvailable" and "NotApplicable" data from a function in which the wind's velocity vector is constructed. The function is

WindFoo[x_, y_] := -WeatherData[{y, x}, "WindSpeed"] {1/Cos[y °], 
   1} Through[{Sin, Cos}[WeatherData[{y, x}, "WindDirection"] °]]

Note: I am not able to display the degree symbol, so it is written out in the above code.

How can the code above be modified to exclude the cases of "NotAvailable" and "NotApplicable" data to work in manner similar to the wind velocity vectors shown over Australia in the following tutorial: http://blog.wolfram.com/2009/01/13/visualizing-weather-patterns-in-mathematica-7/

Answer 1

This uses some version 10 functions:

Getting the data (takes time as I could have set this up better):

dat = Table[{i, j, 
    WeatherData[{j, i}, #] & /@ {"WindSpeed", "WindDirection"}}, {i, 
    113, 153, 4}, {j, -43, -11, 4}];

Processing:

data = Cases[{{#1, #2}, -QuantityMagnitude[#3[[1]]] {1/Cos[#2 Degree], 
        1} Through[{Sin, Cos}[#3[[2]]]]} & @@@ 
    Flatten[dat, 1], {_, {_?NumericQ, _?NumericQ}}];

Visualizing:

ls = ListStreamPlot[data];
arrows = Cases[ls, Arrow[_], Infinity];
GeoGraphics[{Purple, arrows}]

Please look at documentation as must reverse coordinates if you choose to use GeoPosition. Obviously can modify to suit stylistic preferences and other needs.

UPDATE

Just to 'check' this with another source (noting I am not sure what WeatherData uses as default, I just assumed it used the most recent data): From http://www.windfinder.com/weather-maps/forecast/australia#4/-24.69/\ 135.79

or a little before:

These seem somewhat similar to the Mathematica generated map.

Exploiting the hard work of Gleno and appropriately adapting code from the documentation of WindVectorData I believe yields the "correct" wind stream plot consistent with the check above (the issue is transforming the wind velocity components):

windVelocity[
  coords_] := {Reverse@
   coords, -Reverse@QuantityMagnitude[WindVectorData[coords]] {1/
     Cos[coords[[1]] \[Degree]], 1}}

wvd[country_String, date_, n_] := 
  Module[{box, latmin, latmax, longmin, longmax, pts, sample, stream, 
    latstep, longstep, ls, arrows}, 
   box = {Min@#, Max@#} & /@ 
     Transpose@
      Cases[CountryData[country, "Polygon"], {_Real, _}, Infinity];
   {{latmin, latmax}, {longmin, longmax}} = box + {{-2, 2}, {-2, 2}};
   latstep = (latmax - latmin)/(n - 1);
   longstep = (longmax - longmin)/(n - 1);
   pts = Flatten[
     Table[{lat, long}, {lat, latmin, latmax, latstep}, {long, 
       longmin, longmax, longstep}], 1];
   sample = windVelocity /@ pts;
   stream = Cases[sample, u_ /; FreeQ[u, "NotAvailable"]];
   ls = ListStreamPlot[stream];
   arrows = Cases[ls, Arrow[_], Infinity];
   GeoGraphics[arrows, GeoBackground -> "ReliefMap"]
   ];

This leads to:

gr = wvd["Australia", Today, 10]

I guess the answer my friend is blowing in the wind...

Answer 2

NB I think there's a problem with this code, please see this answer by ubpdqn for more information. I'll update these snippets when I get the chance.

Yours is a question with many possible interpretations. I've chosen the interpretation that was most fun for me to play with, so...

ang = 20; (* divide the world into chunks of this size *)
pts = Flatten[Table[{lat, long}, {lat, -80, 80, ang}, {long, -180, 180, ang}], 1];
sample = WindVectorData[Replace[pts, c_ :> {c, Today}, 1]]; (* this takes forever *)

stream = Cases[Transpose@{Reverse/@pts, sample}, {_, Except[_Missing]}] /. 
    Quantity[v_, _] :> v /. Interval[lst_] :> Mean /@ lst;

(* This function crashes kernel very often, but in interpolation, I assume *)
streamplot = ListStreamPlot[stream, StreamColorFunction -> "SolarColors", 
   InterpolationOrder -> None, StreamPoints -> Fine];

mpoly[poly_] := poly /. GeoPosition[pts_] :> pts 
                     /. {lat_Real, long_} :> {long, lat} 
                     /. Polygon[data_] :> Style[Polygon[data], Antialiasing -> True];

world = Graphics[{RGBColor[0.38, 0.77, 1.]}~Join~
    Table[mpoly@CountryData[c, "SchematicPolygon"], {c, CountryData[]}]];
g = Show[world, streamplot, ImageSize -> 1400]

If you have Mathematica 10, the graph g should be:

Puting it together, you can have a function like

PlotWindBox[country_String, date_, n_] := 
  Module[{box, latmin, latmax, longmin, longmax, pts, sample, stream, 
    streamplot, mpoly, target, world, latstep, longstep},

   box = {Min@#, Max@#} & /@ 
     Transpose@
      Cases[CountryData[country, "Polygon"], {_Real, _}, Infinity];
   {{latmin, latmax}, {longmin, longmax}} = box + {{-2, 2}, {-2, 2}};


   latstep = (latmax - latmin)/(n - 1);
   longstep = (longmax - longmin)/(n - 1);
   pts = Flatten[
     Table[{lat, long}, {lat, latmin, latmax, latstep}, {long, 
       longmin, longmax, longstep}], 1];

   sample = WindVectorData[Replace[pts, c_ :> {c, date}, 1]];
   stream = 
    Cases[Transpose@{Reverse /@ pts, sample}, {_, 
        Except[_Missing]}] /. Quantity[v_, _] :> v /. 
     Interval[lst_] :> Mean /@ lst;
   streamplot = 
    ListStreamPlot[stream, StreamColorFunction -> "SolarColors", 
     StreamPoints -> Fine, InterpolationOrder -> None];
   mpoly[poly_] := 
    poly /. 
       GeoPosition[pts_] :> pts /. {lat_Real, long_} :> {long, 
        lat} /. Polygon[data_] :> 
      Style[Polygon[data], Antialiasing -> True];
   target = 
    Graphics[{RGBColor[0.38, 0.77, 1.], 
      mpoly@CountryData[country, "Polygon"]}];
   world = 
    Graphics[{GrayLevel[0.84]}~Join~
      Table[mpoly@CountryData[c, "Polygon"], {c, 
        Cases[CountryData[], Except[Entity["Country", country]]]}]];
   Show[world, target, streamplot, ImageSize -> 1400, 
    PlotRange -> {{longmin, longmax}, {latmin, latmax}}]
   ];

Call it with PlotWindBox["Australia", Today, 10] (10 means that there are 10^2 sample points, try with fewer sample points if you are impatient like me) and obtain:

Answer 3

I feel a bit foolish posting such a simple answer after the two complete examples above, but I think it bears pointing out that version 10 includes a function specifically for removing Missing values:

DeleteMissing: