This could be understood as a physics problem. From mechanics we know that the center of mass is at a point of high symmetry, and if there are symmetry axes (under rotation), then the inertia tensor with respect to the center of mass will have symmetry axes (maybe not all of them - there may be more than there are dimensions) among its principal axes (i.e. eigenvectors).
So given an image, we first convert it to a 2d grid in which each pixel corresponds to a Cartesian coordinate location, and the pixel value is the equivalent of the mass. You have to decide how to translate the color values of the image into a mass distribution. The easiest is to use the brightness.
Given this matrix, the center of mass coordinate com
can be obtained by forming the mass-weighted average of the pixel positions. With that, we can calculate the two-by-two matrix inertia
using the standard formula for the inertia tensor.
The principal axes are obtained from Eigenvectors[inertia]
. Since inertia
is a real symmetric matrix by definition, the eigenvectors are orthogonal. One of the will be the symmetry axis, if there is one.
I'll program the moment of inertia calculation in the most straightforward way without optimizing it (e.g., I don't use the symmetry of inertia
because it's not worth the time). Here is the function that takes a raster image and plots it together with the principal axes relative to the center of mass that were found:
findSym[img_] :=
Module[{imd, nx, ny, points, weights, mass, com, pointsCom, diag,
inertia, v1, v2},
imd = Transpose@
ImageData[RemoveAlphaChannel@ColorConvert[img, "Grayscale"]];
{nx, ny} = Dimensions[imd];
points = Tuples[{Range[nx], Reverse@Range[ny]}];
weights = Flatten[imd];
mass = Total[weights];
com = Total[points weights]/mass // N;
pointsCom = (Subtract[#, com] & /@ points) weights;
diag = Total[(#.#) & /@ pointsCom ];
inertia =
Table[KroneckerDelta[i, j] diag -
Total[(#[[i]] #[[j]]) & /@ pointsCom], {i, 1, 2}, {j, 1, 2}];
{v1, v2} = Eigenvectors[inertia];
Graphics[{Inset[Image[Transpose@imd], Automatic, Automatic,
Scaled[1]], Red, Point[com], Blue,
Line[{{com, nx/3 v1 + com}, {com, nx/3 v2 + com}}]},
PlotRange -> {{1, nx}, {1, ny}}]]
The given input is converted to a grayscale image without alpha channel, so we can extract from it a two-dimensional array of scalar (non-negative) values imd
that represent our mass distribution.
In the output, I reconstruct the image from imd
and add a red dot for the center of mass, and two blue lines showing the principal axes.
Now to test it, I'll import one of your images and rotate it arbitrarily:
g = Import["https://i.stack.imgur.com/lSEfx.gif"];
im = ImageResize[ImageRotate[g, -Pi/3], 150]

Here is the calculation:
findSym[im]

In more complicated pictures, you'd want to pre-process your image by filtering maybe a certain color or selecting only certain morphological components, but that is too wide a field as to be easily covered here.
It would probably be more elegant to reformulate the inertia
definition in terms of Moment
, but I'd do that only if you want to apply this idea to larger images. In general, you probably won't want very high resolution for the analysis anyway because the visual symmetries in art are determined by the coarse features.
Edit
After writing this, I also remembered PrincipalComponents
which finds symmetry axes of an image as well. There is an example in the section "Neat Examples" that does an alignment with the symmetry axes very similar to what the home-made approach above achieves.
The reason the inertia
tensor and PrincipalComponent
approaches are related is that moments of inertia are basically a kind of mass-weighted "variance" of the position relative to an axis. However, if we put a smooth mass-weighting into the input of PrincipalComponents
, we'll get an output that's distorted in shape according to this weighting. Therefore, the example in the documentation uses a binarized input image.
Here is that procedure, applied to the above example:
shape = Position[ImageData[Binarize[im]], 1, {2}];
ListPlot[PrincipalComponents[N@shape], AspectRatio -> Automatic,
PlotRange -> All]

So my inertia
tensor approach above has some advantages because it can work more easily with a smooth "mass" distribution.
Edit 2
As mentioned in the comment, findSym
relies on you to define a mass measure first, and it interprets the image that you pass to it by considering the brighter regions high in mass and the black regions to have zero mass. The latter means the region isn't part of the object for the purposes of the moment of inertia.
To show that this method does not rely on the object's shape as such (i.e., being a rectangle, e.g.), here I mutilate the test image, crop it and throw it onto a black background that bleeds into it, preserving only its internal reflection symmetry along the chandelier:
im = Rasterize[
Graphics[Inset[
ImageRotate[
ImageCrop[
SetAlphaChannel[g,
ColorNegate@
Binarize@
Erosion[Colorize[ClusteringComponents[Binarize@g]],
2]], {226, 270}, Bottom], Pi/3], {Right, Top}, {Right,
Top}], Background -> Black
], "Image"]

findSym[im]

This shows that the symmetry has been identified correctly even though the overall image doesn't have symmetry.
HistogramDistribution
of orientations will be able to locate the center of symmetry. Doesn't this just find the orientations without placing them at the center of symmetry? Maybe I'm overlooking something... for the orientations, it looks like a clever idea though. $\endgroup$