1) The dark image with few groups of high-brightness pixels and some amount of noise around it.
2) Number of clusters.
If the number of clusters is known (like here)
You may use Lloyd's clustering 
The idea is as follows: it optimizes a set of cluster centers $p_i$:
Initialize the p_i's with an initial guess, or randomly For each iteration: Compute the cluster associated with each p_i, (the cluster is the set of points nearer to p_i than to the other p_j's) Move each p_i to the weighted centroid of its cluster
For an image, the iteration can be implented as follows, computing the mass m_i and the centroid g_i of each cluster:
For each i m_i = 0 g_i = (0,0) For each pixel (x,y) of the image let i denote the index of the center p_i nearest to (x,y) m_i = m_i + pixel_intensity(x,y) g_i = g_i + pixel_intensity(x,y) * (x,y) For each i p_i = (1/m_i)*g_i
Since the number of clusters is small, you can find the nearest p_i using a simple loop. If you have a higher number of sites, you may either use a kd-tree, or compute the Voronoi diagram of the sites and iterate on the pixels of each Voronoi cell.
I used this algorithm to cluster the colors of a rubics cube acquired by a lego color sensor, and it works reasonably well while being very easy to implement 
If the number of clusters is unknown then the problem is much more difficult.
You may use "mean shift clustering" , that will apply a filter-like operation to the image, and make the "modes" appear. It acts like the inverse of a smoothing filter.