I want to create a test image which consists of 2d circular shaped objects each having a Gaussian brightness distribution profile.

To make it simple objects should be the same.

I then put the objects at random positions.

In the resulting composed image I would like to prevent overlapping objects.

An alternative would also be that the brighter pixels dominate.

How could that be done?


enter image description here

The single object:

enter image description here

My code:

(*a single object*)
objectData = GaussianMatrix[20, 0];
object = Image[objectData/Max[objectData]];

(*the test image,here 400 pixels square*)
imageData = Array[0 &, {400, 400}];
image = Image[imageData];

(*some random mean positions*)
randomCoords = Array[{RandomReal[], RandomReal[]} &, 20];

(*image composition*)
ImageCompose[image, object, randomCoords*400]
  • 1
    $\begingroup$ closely related: Adding and summing overlapping Airy (or Gaussian) disks on a plane $\endgroup$ – Kuba Feb 18 '16 at 17:44
  • 3
    $\begingroup$ GaussianFilter[RandomImage[BernoulliDistribution[.0002], {300, 300}], 10] // ImageAdjust $\endgroup$ – Dr. belisarius Feb 18 '16 at 17:46
  • 2
    $\begingroup$ Strictly speaking, since Gaussians don't have compact support, the objects will always be overlapping unless you truncate the Gaussians. It might become clearer what the best answer would be if you specify what the desired purpose is. $\endgroup$ – Jens Feb 18 '16 at 18:51
  • $\begingroup$ The objects should represent simulated spherical micro-particles that are illuminated with a laser plane and recorded perpendicular to the laser plane with a camera. The scattered light from the particles is in the ideal case Gaussian distributed. $\endgroup$ – mrz Feb 19 '16 at 14:16
  • $\begingroup$ The code shown appears to be riddled with name changes e.g. image/img, object/particle. $\endgroup$ – Daniel Lichtblau Feb 19 '16 at 17:36


Belisarius suggested in comments a method using image processing (i.e. GaussianFilter) to generate the image. I'd modify his approach to use RandomChoice to generate the initial spots, which seems more direct to me:

img = Image@RandomChoice[{10000, 1} -> {0, 1}, {400, 400}];
ImageAdjust@GaussianFilter[img, 15]

Mathematica graphics

I am going to create your object as the PDF of a binormal distribution with means equal to the $(x,y)$ position of the object, and equal standard deviations (to give a circular particle image) whose values (8) have been chosen to be pleasing with an image size of 400 units (of course, adjust to taste).

Clear[object, positions, imagefunction]

object[px_, py_] = PDF[BinormalDistribution[{px, py}, {8, 8}, 0], {x, y}];

Now generate some random positions, apply the object function to generate those objects, then sum all those distribution PDF to obtain a parametric description of a 3D surface:

positions = RandomInteger[{1, 400}, {20, 2}];
imagefunction = Plus @@ object @@@ positions;

Finally use DensityPlot to slice that surface and obtain your image:

DensityPlot[imagefunction, {x, 1, 400}, {y, 1, 400},
  PlotPoints -> 150, ColorFunction -> GrayLevel,
  Frame -> None, PlotRange -> All

Mathematica graphics

To get an actual image, rather than a Graphics object, use e.g. Rasterize, Export ...:

| improve this answer | |
  • $\begingroup$ Thank you very much for your help. Also for your solution I have the same problem that the intensities of overlapping objects are higher than for a single one. What I would like to have: Intensity(pixel) = Max [ object1Intensity(pixel), object2Intensity(pixel) ]? (see comment to @Dr belisarius) $\endgroup$ – mrz Feb 19 '16 at 13:55

One way to proceed is to make an image that contains a single white pixel at each desired center point and then to convolve the small particle with the image:

particle = Import["http://i.stack.imgur.com/gi6U1.png"];
randomCoords = Array[{RandomReal[], RandomReal[]} &, 20];
img = Image[Array[0 &, {400, 400}]];
points = ReplacePixelValue[img, randomCoords*400 -> 1];
ImageConvolve[points, particle]

enter image description here

| improve this answer | |

ImageAdjustwon't work well when there are superpositions,

n = 400;
gf[i_Image] := GaussianFilter[i, {n/20, n/80}]
ImageAdjust [ ri= gf@RandomImage[BernoulliDistribution[50/n^2], {n, n}]]

Mathematica graphics

Use Clip/Rescale instead:

resc = Max@ImageData@gf@Image@SparseArray[{n/2, n/2} -> 1, {n, n}];

Image@Clip@Rescale[ImageData@ri, {0, resc}]
(* or, mostly the same but slower*)
ImageApply[Clip@Rescale[#, {0, resc}] &, ri]

Mathematica graphics

| improve this answer | |
  • $\begingroup$ Great solution, only one question: why is the brightness of the overlapping objects higher than the rest? $\endgroup$ – mrz Feb 18 '16 at 20:17
  • $\begingroup$ @mrz In the rescaled code it isn't "higher", what happens is that the maxima are extended (in the non.rescaled code it is, that's why I suggested to rescale) $\endgroup$ – Dr. belisarius Feb 18 '16 at 20:34
  • $\begingroup$ When I measure the maxium brightness or histogram of overlapping objects, then their brightness is higher than of non-overlapping objects goo.gl/jIy6ND. $\endgroup$ – mrz Feb 19 '16 at 11:31
  • $\begingroup$ Oh, I see,you're using ImageData. Let me check $\endgroup$ – Dr. belisarius Feb 19 '16 at 12:08
  • 1
    $\begingroup$ @mrz Yes, But it doesn't look nice to me is = ImageAdjust /@ gf /@ Image /@ (SparseArray[{#} -> 1, {n, n}] & /@ RandomInteger[{1, n}, {50, 2}]); ImageApply[Max, is] $\endgroup$ – Dr. belisarius Feb 19 '16 at 13:28

Here is a pedestrian way to avoid overlaps. I'll show it "by hand" although it could be automated. I believe this or something similar came up once before on this forum but I cannot seem to locate it.

The original setup is as follows.

rad = 20;
scale = 400;
objectData = GaussianMatrix[rad, 0];
object = Image[objectData/Max[objectData]];
imageData = Array[0 &, {scale, scale}];
image = Image[imageData];
randomCoords = RandomReal[1, {20, 2}];
image1 = ImageCompose[image, object, randomCoords*scale]

enter image description here

To get something similar but without overlaps we can iteratively use Nearest. I'll show one approach although there are others. The idea is to select only points that have no neighbor too close. This is inefficient insofar as we could take one from any pair that overlap, but implementing that means in effect changing the NearestFunction in a loop.

randomCoords = Array[{RandomReal[], RandomReal[]} &, 40];
nf = Nearest[randomCoords];
randomCoords2 = 
   With[{nbr = nf[#, {2, 2/Sqrt[scale]}]}, Length[nbr] == 1] &];
randomCoords3 = Array[{RandomReal[], RandomReal[]} &, 30];
nf2 = Nearest[Join[randomCoords2, randomCoords3]];
randomCoords4 = 
   With[{nbr = nf2[#, {2, 2/Sqrt[scale]}]}, Length[nbr] == 1] &];

(* Out[246]= 11

Out[250]= 9 *)

We again have 20 points as luck would have it. So let's look at the placement this time.

image2 = ImageCompose[image, object, 
  Join[randomCoords2, randomCoords4]*400]

enter image description here

| improve this answer | |

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.