Calculate the size of the moon
Mathematica 126 119 109 bytes
Mathematica can measure the elongation of a component in an image.  A full moon, being perfectly symmetric, has an elongation of 0, on a scale of 0 to 1.
A diminishing moon becomes increasingly elongated, to a maximum of roughly 0.8.
0.998 -0.788 x-0.578 x^2 was the empirically determined model (based on the large photos) for `predicting the fullness of the moon (by area), given its elongation.
I adjusted the model to 1- 0.788 x -0.578 x^2 so that with exactly zero elongation (full moon) the model will return 1 for the pixel scale factor.  It saves 4 bytes and still stays within the accuracy limits.
This model is used for any size images.  The moon image does not need to be centered. It also does not need to cover a fixed proportion of the photo.
Here are the data points (elongation, displayedMoonPixels/fullMoonPixels) for the large images and the parabolic model that was generated to fit the data. Linear models fit ok, but the quadratic model is dead on, within limits (see below).
Here the data are from the large pictures.  So is the model


Below, the data (the red points) are from the small pictures.  The model (the blue curve) is the one generated by the large pictures, the same one as displayed above.
The smallest crescent has 7.5% the area of a full moon.  (The smallest crescent among the large photos is 19% of a full moon.)  If the quadratic model had been based on the small photos the fit below would be better, only because it accommodated the small crescent.  A robust model, that would stand up under a wide range of conditions, including very small crescents, would be better made from a larger variety of pictures.
The closeness of fit shows that the model was not hard-coded for the given pictures. We can be fairly certain that the elongation of a moon is independent of the size of the photo, as one would expect.

f takes the image, i, as input and outputs the predicted size of the full moon, in pixels.  It works for off-center shots.
As the data below show, it all of the test cases except one.
The moons were arranged from full to most diminished.
i_~c~t_ := Max@ComponentMeasurements[i, t][[All, 2]];
f@i_ := i~c~&quot;Count&quot;/(1 - 0.788 x - 0.578 x^2 /. x -&gt; i~c~&quot;Elongation&quot;)

More than one image component may turn up in a photo.  Even a single pixel separated from the others will be considered a distinct component.  For this reason, it is necessary to search &quot;all&quot; components, to find the one that has the greater number of pixels.  (One of the small photos has more than one image component.)
Large pictures
Predictions of moon size made from the large photos were uniformly accurate.
{&quot;predicted size of full moon&quot;, f[#] &amp; /@ large}
{&quot;accuracy&quot;, %[[2]]/16416}


{&quot;predicted sizes of full moon&quot;, {16422., 16270.9, 16420.6, 16585.5, 16126.5, 16151.6}}
{&quot;accuracy&quot;, {1.00037, 0.991161, 1.00028, 1.01033, 0.982367, 0.983891}}


Small pictures
Predictions of moon size made from the small photos were uniformly, with one great exception, the final picture.  I suspect the issue stems from the fact that the crescent is very narrow.
{&quot;predicted sizes of full moon&quot;, f[#] &amp; /@ small}
{&quot;accuracy&quot;, %[[2]]/10241}


{&quot;predicted sizes of full moon&quot;,{10247.3, 10161., 10265.6, 10391., 10058.9, 7045.91}}
{&quot;accuracy&quot;, {1.00061, 0.992192, 1.0024, 1.01465, 0.982221, 0.68801}}