# 3d points projected to y-z plane and colorized according x-velocity using a grid for averaging

This Question is related to this solved by @kglr and also this solved by @Henrik Schumacher.

I have a set of 3d positions and velocities of particles, and want to colorize the projected positions to the xz plane according to the x-velocity.

kglr solved this by usinge the following code:

SeedRandom[1];
posdata = RandomReal[10, {1000, 3}];
veldata = RandomReal[{-0.5, +0.5}, {1000, 3}];

clfun[val_] := Blend[{Blue, Green, Yellow, Red}, val];

xvelzypos = posdata;
xvelzypos[[All, 1]] = veldata[[All, 1]];

{xmin, xmax} = MinMax@veldata[[All, 1]];

legend = BarLegend[{clfun[Rescale[#, {xmin, xmax}]] &, {xmin, xmax}},
LegendLabel -> "vx", LegendMarkerSize -> {20, 300},
LabelStyle -> {FontFamily -> "Calibri", FontSize -> 15}];

xvzyplot =
Show[Graphics[{PointSize[0.03],
clfun[Rescale[#[[1]], {xmin, xmax}]],
Point[{#[[3]], #[[2]]}]} & /@ xvelzypos, Frame -> True,
FrameLabel -> {{"y", ""}, {"z", ""}}, AspectRatio -> Automatic,
PlotRange -> All,
BaseStyle -> {FontWeight -> "Bold", FontSize -> 20},
ImageSize -> 400]];

Legended[Graphics[
Inset[xvzyplot, Scaled[{.5, .5}], Automatic, Scaled[1]],
AspectRatio -> ImageAspectRatio@xvzyplot], legend]


The result is:

My question:

I would like to use a certain grid on the yz plane with values dy and dz. Each grid element should be color coded with the mean velocity of all particles that are projected to the corresponding grid element.

I tried to use parts of the solution of @Henrik Schumacher, but did not succeed.

How can be my problem solved?

Here is one a out of probably a million possibilities to achieve that:

SeedRandom[1];
posdata = RandomReal[10, {1000, 3}];
veldata = RandomReal[{-0.5, +0.5}, {1000, 3}];
ymin = 0.;
ymax = 10.;
zmin = 0.;
zmax = 10.;
dy = .5;
dz = .25;
m = Quotient[ymax - ymin, dy];
n = Quotient[zmax - zmin, dz];
boxcenters = Tuples[{
Range[ymin + 0.5 dy, ymax - 0.5 dy, dy],
Range[zmin + 0.5 dz, zmax - 0.5 dz, dz]
}];


If the aspect ration of the boxes is not equal to 1, the trick is to rescale the positions so that we have quadratic boxes. This way, we can employ Nearest with DistanceFunction -> ChessboardDistance to find the correct box for each point.

normalizedboxcenters = Tuples[{Range[.5, m - 0.5], Range[.5, n - 0.5]}];
idx = Flatten@Nearest[
normalizedboxcenters -> Automatic,
Transpose[{(#2 - ymin)/dy, (#3 - zmin)/dz} & @@ Transpose[posdata]]
,
{1, 0.5},
DistanceFunction -> ChessboardDistance
];


Have a index of the box for each point, we build a sparse matrix that will do the summing up by matrix-vector multiplication. Total[A, {2}] tells us, how many points are in each box which also helps us to perform the averaging.

A = SparseArray[
Transpose[{idx, Range[Length[idx]]}] -> 1.,
{Length[boxcenters], Length[posdata]},
0.
];
meanvelocities = (A.veldata)/Clip[Total[A, {2}], {1., \[Infinity]}];
meanxvelocities = meanvelocities[[All, 1]];


Assigning some colors and plotting the result.

cols = ColorData["DarkRainbow"] /@ Rescale[meanxvelocities];
Graphics[
Transpose[{
cols,
Rectangle @@@ Transpose[{
Transpose[Transpose[boxcenters] - 0.5 {dy, dz}],
Transpose[Transpose[boxcenters] + 0.5 {dy, dz}]
}]
}]
]


Of course, you are free to prettify the output and to add labels and stuff. But I guess, you get the idea.

• Thank you so much for this solution ,,, it helps me a lot – mrz Jul 13 '18 at 16:32
• You're welcome. – Henrik Schumacher Jul 13 '18 at 16:33
• I tried your solution with my real data set and have a problem with it. Can you please have short look on it? Thank you. drive.google.com/open?id=1TXoUllMcws-ax1SMQFAd9QMeE4t5jrGv – mrz Jul 23 '18 at 13:39
• Okay, I found a typo of mine. Please try again with boxcenters = Tuples[{Range[ymin + 0.5 dy, ymax - 0.5 dy, dy],Range[zmin + 0.5 dz, zmax - 0.5 dz, dz]}];. – Henrik Schumacher Jul 23 '18 at 13:53
• Not also that you it is possible that you get your image reflected with respect to the main diagonal. If you have pixels that are squares, it might be more efficient with ArrayPlot[ArrayReshape[Rescale[meanxvelocities], {m, n}], ColorFunction -> ColorData["DarkRainbow"]]. – Henrik Schumacher Jul 23 '18 at 13:55

### SparseArray and "TreatRepeatedEntries"

Related Q/As :

vvalues = xvelzypos[[All, 1]];
epsilon = 1*^-10;
indexes = 1 + Floor[(1 - epsilon) 20 Rescale[xvelzypos[[All, {2, 3}]]]];
SystemSetSystemOptions["SparseArrayOptions" -> {"TreatRepeatedEntries"->(Mean[{##}]&)}];
binmeansV = SparseArray[indexes -> vvalues];
SystemSetSystemOptions["SparseArrayOptions" -> {"TreatRepeatedEntries" -> First}];

MatrixPlot[binmeansV, PlotRange -> {{0, 10}, {0, 10}},
BaseStyle -> Opacity[.7], ColorFunction -> clfun,
ImageSize -> 1 -> 40,
FrameTicks -> {{#, #} & @ Transpose[{#, # - .5} & @Range[.5, 10.5, 1]],
{#, #} & @ Transpose[{#, # - .5} &@Range[.5, 10.5, 1]]}]


### BinLists

Using 5 bins each of size 2 for the horizontal and vertical coordinates and a single bin for velocity:

xbinspec = ybinspec = {0, 11, 2};
vbinspec ={-1, 1, 2};
binmeans = Flatten[Map[Mean[Last /@ #] &,
BinLists[RotateLeft /@ xvelzypos, xbinspec , ybinspec , vbinspec], {-3}], {3}][[1]];
matrixplot = MatrixPlot[Riffle[#, #] &[Riffle[#, #] & /@ binmeans],
PlotRange -> {{0, 10}, {0, 10}}, BaseStyle -> Opacity[.7],
ColorFunction -> clfun, ImageSize -> 1 -> 40, DataReversed -> True,
FrameTicks -> {{#, #} &@Transpose[{#, # - .5} &@Range[.5, 10.5, 1]],
{#, #} &@Transpose[{#, # - .5} &@Range[.5, 10.5, 1]]}]


Add the option Epilog -> xvzyplot[[1]] to get

• Your explanation and solution is excellent … I will go through it. My real data consist of about 10^6 to 10^7 data points. Thank you very much for your help. – mrz Jul 14 '18 at 8:22
• I tried your solution with my real data set and have a problem with it. Can you please have look on the following notebook? Thank you. drive.google.com/open?id=1JQMop_lcxI2pL5tUvH9A17KOF3tTV0s5 My data range is not exactly squared, this is one problem. I know that I did not properly set the PlotRange. I tried to adapt it to MinMax of posdata[[All,2]] and posdata[[All,3]] but it did not help. How can I increase the number of row and columns? – mrz Jul 24 '18 at 6:35
• @Lenoil, for the SparseArray method, you can use indexes = 1 + Floor[(1 - epsilon) {numberofxbins, numberofybins} # & /@ Rescale[xvelzypos[[All, {2, 3}]]]] ti have different number of bins for the first and second dimension. – kglr Jul 24 '18 at 7:42
• @L, Btw, you also need to change the PlotRange option from MatrixPlot. – kglr Jul 24 '18 at 7:49