CODEWITHSUNDEEP
pythonnumpy
XaBla
XaBla

Reputation: 127

How to count how many particles in a 2D Gaussian distribution

I've just started using Python, so I'm sorry if I ask for trivial things.

I generated two random Gaussian distributions, and I used them to generate a 2D Gaussian distribution. What I'd like to do now is to plot a graph that represents the number of elements within a circumference of the 2D Gaussian distribution, varying the radius of the circumference (reducing it at each step).

You would be so kind to help me solve the problem. Thank you for taking my post into consideration.

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns

mu1, sigma1 = 0, 0.1                                  # mean and standard deviation
s1 = np.random.normal(mu1, sigma1, 10000)              # generate N randoomly Gaussian points
mu2, sigma2 = 0.8, 0.3
s2 = np.random.normal(mu2, sigma2, 10000)


#ISTOGRAMMA DI DUE DISTRIBUZIONI GAUSSIANE CON DIFFERENTI SET
plt.figure(1)
plt.title('Histogram of a 2D-Gaussian Distribution')
bins1 = plt.hist(s1, 100)
bins2 = plt.hist(s2, 100)
plt.show()

#DISTRIBUZIONE GAUSSIANA 2D
plt.figure(2)
plt.title('2D-Gaussian Distribution')
bins = plt.hist2d(s1, s2, 100)
cb = plt.colorbar()
cb.set_label('counts in bin')
plt.show()

Upvotes: 0

Views: 370

Answers (1)

HitLuca
HitLuca

Reputation: 1280

This should do the trick, assuming that s1 and s2 are the coordinates of some 2D points. If not, the code can be easily changed to match your problem.

First you center the two distributions by subtracting their mean, then you check which of their elements (absolute value) is inside the radius of your circle. You then take a logical and to make sure only to take the points that have both coordinates inside the circle.

radius = 0.1
valid_indexes = np.logical_and(abs(s1 -mu1)<= radius, abs(s2 - mu2) <= radius)

s1_valid = s1[valid_indexes]
s2_valid = s2[valid_indexes]

You have now obtained the points in the distributions that are inside a circle with given radius centered in (mu1, mu2).

[Edit]

As you want to count the number of elements, and not extract them, you can easily do

radius = 0.1
sum(np.logical_and(abs(s1 -mu1)< radius, abs(s2 - mu2) < radius))

[Edit 2]

This plots the number of points for every radius of the circle, starting from limit and reducing it by step until 0

step = 0.025
limit = 1

s1_ca = abs(s1-mu1)
s2_ca = abs(s2-mu2)

points_in_radius = []

radius_values = np.round(np.arange(0, limit, step), 3)[::-1]

for radius in radius_values:
    points_in_radius.append(sum(np.logical_and(s1_ca < radius, s2_ca < radius)))

plt.plot(points_in_radius)
plt.xticks(range(len(points_in_radius)), radius_values, rotation=90)
plt.show()

first i center the distributions and take their abs value. then I create the range of radiuses to use and finally I cycle through them and add the result using my above formula.

This is not the most efficient way to do it, but it works.

Upvotes: 1

Related Questions