Functions

Set of rotation matrices in 3D

rotX[α_]:={{1,0,0},{0,Cos[α],-Sin[α]},{0,Sin[α],Cos[α]}};
rotY[α_]:={{Cos[α],0,Sin[α]},{0,1,0},{-Sin[α],0,Cos[α]}};
rotZ[α_]:={{Cos[α],-Sin[α],0},{Sin[α],Cos[α],0},{0,0,1}};

Funtion to rotate the objects within a graph. 
The input is similar to that of ParametricPlot3D[___] including the Options[___]. 
There are 6 additional arguments:
	{x,y,z} a vector that shifts the object and
	{α,β,γ} three Euler-angles (in the z-y-z-system) that rotate the object

shiftAndRotate[fun_,{θ_,θs_,θe_},{φ_,φs_,φe_},{x_,y_,z_},{α_,β_,γ_},opts___]:=Module[
{
graph,points,rotationMatrix=rotZ[γ].rotY[β].rotZ[α],rotationMatrixT,
normPos,normals,normalVecs
},
rotationMatrixT=Transpose[rotationMatrix];
graph=ParametricPlot3D[ fun,{θ,θs,θe},{φ,φs,φe},opts];
points=graph[[1,1]];
points={x,y,z}+rotationMatrix.#&/@points;
graph[[1,1]]=points;
normPos=Flatten[Position[graph,VertexNormals]];
If[Length[normPos]>0,
normPos=Drop[normPos,-1];
normals=graph[[Sequence@@normPos]];
normalVecs=normals[[-1]];
normalVecs=rotationMatrix.#&/@normalVecs;
normals[[-1]]=normalVecs;
graph[[Sequence@@normPos]]=normals;
];
Return[graph]
]

As a simple example look at this quadrupole moment

shiftAndRotate[ {Sin[φ] Sin[θ], Cos[φ]Sin[θ], Cos[θ]}Abs[SphericalHarmonicY[2,0,θ,φ]], {θ, 0, π}, {φ, 0,2 π},{1,0,0},{10 °,20 °,30 °},MaxRecursion->1,PlotStyle->{Specularity[White,25]},
Axes->True,Mesh->None,
PlotPoints->{25,25},
ColorFunction->Function[{x,y,z,θ,φ},
If[SphericalHarmonicY[2,0,π θ,2π φ]>0,
	If[Sin[2π 6φ]>0,ColorData["SunsetColors"][((Cos[45°])^2)],RGBColor[1,1,1]],
If[Sin[2π 6φ]>0,ColorData["SunsetColors"][((Cos[45°])^2)],RGBColor[0,0,0]]]],
Lighting->{{"Directional",GrayLevel[1],ImageScaled[{1,1,1}]}},NormalsFunction->None,PlotRange->All
]

Arrays

The graphs maybe combined like this

In general one has a list of moment coordinates containing positions and rotations

positionList=Flatten[Table[{{1.5i,1.5j,0},{0,(π/2.),Random[Real,{0,2π}]}},{i,4},{j,4}],1];

graphTable={};Do[
AppendTo[graphTable,shiftAndRotate[ {Sin[φ] Sin[θ], Cos[φ]Sin[θ], Cos[θ]}Abs[SphericalHarmonicY[1,0,θ,φ]], {θ, 0, π}, {φ, 0,2 π},positionList[[n,1]],positionList[[n,2]],
MaxRecursion->2,PlotStyle->{Specularity[White,25]},
Axes->True,Mesh->None,
PlotPoints->{15,15},
ColorFunction->Function[{x,y,z,θ,φ},
If[SphericalHarmonicY[1,0,π θ,2π φ]>0,RGBColor[.8,0,0],RGBColor[0,.6,0] ]],NormalsFunction->Automatic,PlotRange->All
]
],{n,1,Length[positionList]}]

Show[graphTable,PlotRange->All,Lighting->{
{"Point",GrayLevel[1],ImageScaled[{-5,-5,-5}]},
{"Directional",GrayLevel[.5],ImageScaled[{0,1,1}]}
}]

It is also possible to make an orientation dependent coloring (in this case with respect to |z|)

positionList2=Flatten[Table[{{1.5i,1.5j,0},{0,(π/16)(i+j),(π/11)(i-j)}},{i,0,8},{j,0,8}],1];

graphTable2={};Do[
AppendTo[graphTable2,shiftAndRotate[ {Sin[φ] Sin[θ], Cos[φ]Sin[θ], Cos[θ]}Abs[SphericalHarmonicY[2,0,θ,φ]], {θ, 0, π}, {φ, 0,2 π},positionList2[[n,1]],positionList2[[n,2]],
PlotStyle->{Specularity[White,25]},
Axes->True,Mesh->None,
PlotPoints->{15,15},
MaxRecursion->0,
ColorFunction->Function[{x,y,z,θ,φ},
If[SphericalHarmonicY[2,0,π θ,2π φ]>0,
	If[Sin[2π 3φ]>0,ColorData["Rainbow"][((Cos[positionList2[[n,2,2]]])^2)],RGBColor[1,1,1]],
If[Sin[2π 3φ]>0,ColorData["Rainbow"][((Cos[positionList2[[n,2,2]]])^2)],RGBColor[0,0,0]]]],
NormalsFunction->Automatic
]
],{n,1,Length[positionList2]}]

Show[graphTable2,PlotRange->All,Lighting->{
{"Directional",GrayLevel[.6],ImageScaled[{1,1,1}]},
{"Directional",GrayLevel[.6],ImageScaled[{0,1,1}]}
}]