Modified ListPlot, how to speed it up

Question

I use Mathematica almost every day for several jobs (computations, graphics, etc.). I admit that this software contains much more options (capabilities) than one user may need. Well, not always since I hit a snag. Let me be more specific:

I have a data set which contains hundreds of thousands of sub-lists composing of four items; the first two are the $(x,y)$ coordinates the third is a real number $t$ and the four is always an integer $i$. For the purpose of this question let's generate a random data

n = 50000;
data = Table[{RandomReal[], RandomReal[], RandomReal[{10^-2, 10^4}], 
RandomInteger[{0, 2}]}, {i, 1, n}];

OK, now we want to create a two-dimensional plot depicting all the $(x,y)$ points with $i \ne 0$ and assign to each point a color according to the value $log_{10}f$. One may reasonably argue, that for what I describe here I need either ListDensityPlot or ListContourPlot. Unfortunately, both modules come by default with a nasty drawback which is fatal in my case. The problem is that both modules spread the colors or in other words, they try to join the colors between the points. There should be an option like Joined -> False but it isn't. Perhaps this should be forwarded to Wolfram, so as to be considered for the new version of the program.

So, for creating my plot I have to push the program to the limits.

First, we filter the data, we select only those sub lists with $i ne 0$ and we define a color function (for those which $n = 0$ we assign gray color).

valrange = {-2, 4};
data[[All, 3]] = Rescale[Log10[data[[All, 3]]] // N, valrange];
colfunc[x_, cf_] := If[x[[4]] == 0, Gray, ColorData[cf][1 - x[[3]]]];

Then, more labor work is needed to create the corresponding pallet explaining the colors

Clear[colorbar]
colorbar[{min_, max_}, colorFunction_: Automatic, divs_: 150] := 
DensityPlot[y, {x, 0, 0.1}, {y, min, max}, AspectRatio -> 10, 
PlotRangePadding -> 0, PlotPoints -> {2, divs}, MaxRecursion -> 0, 
Frame -> True, 
FrameLabel -> {{None, Row[{Subscript["log", "10"], "(", Subscript["t", "esc"], ")"}]}, {None, None}}, 
LabelStyle -> Directive[FontFamily -> "Helvetica", 17], 
FrameTicks -> {{None, All}, {None, None}}, 
FrameTicksStyle -> Directive[FontFamily -> "Helvetica", 15, Plain], 
ColorFunction -> colorFunction]

and finally use ListPlot to do the job

With[{opts = {ImageSize -> {Automatic, 550}}, cf = "Rainbow"}, 
 Row[{Show[
ListPlot[List /@ data[[All, {1, 2}]], 
 PlotStyle -> ({PointSize[0.005], colfunc[#, cf]} & /@ 1. data), 
 AspectRatio -> 1, Frame -> True, RotateLabel -> False, 
 Axes -> None, FrameTicks -> True, FrameLabel -> {"x", "y"}, 
 LabelStyle -> Directive[FontFamily -> "Helvetica", 20], 
 ImagePadding -> {{60, 20}, {60, 20}}, opts, 
 PlotLabel -> 
  Style["Test 2", FontSize -> 20, Black, 
   FontFamily -> "Helvetica"]], PlotRange -> 0.9, 
PlotRangeClipping -> True], 
Show[colorbar[valrange, ColorData[cf][1 - #] &], 
ImagePadding -> {{20, 60}, {60, 40}}, opts]}]]

This piece of code works nice when the number of points is relatively small, let's say $n < 10^4$. My real data however, correspond to $n > 5 \times 10^4$. My CPU (Dual Core at 2.2GHz) needs about 15 minutes for this task. On the other hand, when I enter the same amount of data in already built-in modules (i.e use ListPlot to plot only the coordinates) the required time does not exceed 60 secs.

So, the question: how can I speed up this procedure and lower the CPU time down to let's say no more than 2 minutes? (For testing keep $n = 50000$ which corresponds roughly to a real data set).

Answer 1

This quick test of @ssch's suggestion in the comments:

n = 100000;
data = Table[{RandomReal[], RandomReal[], RandomReal[{10^-2, 10^4}], 
    RandomInteger[{0, 2}]}, {i, 1, n}];
colfunc[x_, cf_] := If[x[[4]] == 0, Gray, ColorData[cf][1 - x[[3]]]];

Graphics[{PointSize[Tiny], 
  Point[#[[1 ;; 2]], VertexColors -> colfunc[#, "Rainbow"]] & /@ data}] 

takes about 10 seconds on my 2011 iMac. For no good reason:

(I think I've discovered dark matter...)

Answer 2

I would pre-process the data:

n = 50000;
data = Table[{RandomReal[], RandomReal[], RandomReal[{10^-2, 10^4}], RandomInteger[2]}, 
             {i, 1, n}];

data[[ ;; , 3]] = Map[Blend["Rainbow", 1 - #] &,
                      Rescale[Log10[data[[ ;; , 3]]], {-2, 4}, {0, 1}]
                     ];

data[[ ;; , 3]] = ReplacePart[data[[ ;; , 3]], Position[data[[ ;; , 4]], 0] -> Gray];

bar = ArrayPlot[Array[# &, {121, 1}], Frame -> {{None, True}, {None, None}}, 
         FrameTicks -> {{None, Transpose[{#, Reverse@Rescale[#, {1, 121}, {-2, 4}]
                        } &@Range[1, 121, 20]]}, {None, None}}, 
         AspectRatio -> 10, ColorFunction -> "Rainbow", ImageSize -> {Automatic, 400}];


Row[{ Graphics[Point[data[[ ;; , ;; 2]], VertexColors -> data[[ ;; , 3]]],
               BaseStyle -> AbsolutePointSize[5], ImageSize -> 400],
     bar}]