moritzhuetten / dmbounds

A database for indirect Dark Matter limits
Other
5 stars 1 forks source link

Casting of limits #21

Closed micheledoro closed 11 months ago

micheledoro commented 1 year ago

Hi @kumark0 I propose we continue here the discussion about the casting of limits. As soon as you have updates, please post them here for proper book-keeping.

kumark0 commented 1 year ago

Dear, @micheledoro sir,

I have followed the steps as you have mentioned. I took the limits for bb channel from the paper and multiplied them with ratio of dN/dx(bb/tautau). After that i got a series of values of tautau limit at a given dark matter mass for all energy fraction x=2E/m and then i took the maximum of it to take it as the limit of tautau channel at the given dark matter mass. I did the same in a loop for all the dark matter masses whose bb limits have been given to calculate the corresponding limit for tautau.

kumark0 commented 1 year ago

[HDM limits - Jupyter Notebook.pdf](https://github.com/moritzhuetten/dmbounds/files/11768501/HDM.limits.-.Jupyter.Notebook.pdf)

micheledoro commented 1 year ago

Considering the ratio between the two spectra is around 1, there is not way the two compared can be so far off. There must be some bug in the notebook. If you attach the .ipynb or git pull it I can have a look

kumark0 commented 1 year ago

Ok, sir. I have attached the kink for the notebook file.

https://github.com/kumark0/HDM_limits/blob/main/HDM%20limits.ipynb

micheledoro commented 1 year ago

Hi @kumark0. I had a look at the notebook, it appears to me 1/ because we need to compare with published data you need both the 'bb_published' and 'tautau_published' and then you compare 'tauta_guessed' with 'tautau_published' 2/ considering the ratio is between 0.1 and 10 there must be something wrong with the comparison plot. From my quick view it looks like you are mixing masses in TeV and GeV but I am not sure this is the case 3/ you can create a folder notebook in this repository and 'git push' your notebooks there if you prefer. It would make the code sharing easier. 4/ I can run HDMSpectra only from a specific folder, do you know how to be able to import it from any location?

kumark0 commented 1 year ago

Dear @micheledoro sir,

(1,2) I'll check the masses scales once again and do the correction as you have mentioned and inform you.

(4) Once you install it using "python setup.py install" in the HDMSpectra folder you have downloaded, it can be accessed anywhere from the folder.

(3) Also, I'll create the repo folder as you have mentioned.

kumark0 commented 1 year ago

Dear @micheledoro sir, I created the notebook folder as u have mentioned in the last text and added the HDM_limit ipynb file there. I checked the masses scales, and everything seems to be in GeV Scale. Can you please check the updated notebook file?

micheledoro commented 1 year ago

Hi @kumark0 can you write a written report here with images and comments so that I can comment inline? I believe that just going through notebook is not the most efficient solution. I would like to see:

Images can be easily copy-paste in this markdown messages. Thanks.

kumark0 commented 1 year ago

sure sir, I'll make a slide for the same having each steps i did and share it with u.

micheledoro commented 1 year ago

There is no need of slides, just copy the relevant images here in order to keep all here without need of external files

kumark0 commented 1 year ago

ok sir

kumark0 commented 1 year ago

Dear @micheledoro sir,

(1)

Original bb and tau tau limits have been taken from the GitHub repository you have shared and referred to the paper "Optimized dark matter searches in deep observations of Segue 1 with MAGIC".

Link: https://github.com/moritzhuetten/dmbounds/blob/main/dmbounds/bounds/magic/magic_2014_segue1_dec_bb.ecsv

mDM =np.array([1.105295e3 ,1.365418e3,1.783231e3 ,2.517480e3 ,3.757329e3,4.961947e3,7.162629e3 ,9.889380e3, 13.96135e3 ,19.93044e3]) bb_limits=np.array([1.736478e24,2.516671e24,3.544759e24,5.596723e24,8.836514e24,1.110336e25,1.435571e25,1.803841e25, 2.140811e25,2.399746e25])

tau-tau limits taken from : https://github.com/moritzhuetten/dmbounds/blob/main/dmbounds/bounds/magic/magic_2014_segue1_dec_tautau.ecsv

Similarly sigmaV limit were taken from the ref :

https://github.com/moritzhuetten/dmbounds/blob/main/dmbounds/bounds/magic/magic_2014_segue1_ann_bb.ecsv

https://github.com/moritzhuetten/dmbounds/blob/main/dmbounds/bounds/magic/magic_2014_segue1_ann_tautau.ecsv

###############################################################################################

(2)

I calculated the ratio of bb/tau-tau for each DM mass in for loop for all the energy fraction values. Hence the output ratio is a numpy array for all the energy fractions x =2E/m. Here is the plot for same :

133

Then, I multiplied the bb limit by the ratio obtained and took its maximum value. Initially, I have a ratio for all energy fraction values. After multiplying to bb limit, I got an array of limit values for tautau for all energy fractions and I took the maximum value of it for the lifetime limit of tau tau and the minimum value of it for the case of limit of annihilation. Here is the plot showing the distribution of the lifetime limit of tau tau i obtained after the multiplication of ratio with bb limit.

image

The plot shows the distribution of tautau limit for each of dark matter masses.

###########################################################################################

(3)

Once i had them, i plotted them with the published limits and here are the results :

image

Same process i did for the case of annihilation, and here are the plots for ratios, distribution of limits over dark matter masses, and final plots for comparing the calculated and published limits.

image

image

Discussion:

i had a sequence of values of ratio for all energy fractions, and when I multiplied them with a limit of bb, I got a sequence of values for tautau limit for each of the dark matter masses. I took the maximum of it as a limit for the decay case and the minium of it as a limit for the annihilation case after multiplying the corresponding bb decay and annihilation limits. Can you please check that part?

micheledoro commented 1 year ago

Hi @kumark0 apologies for the late reply. Comments

Annihilation/decay

There is no need to use both decay and annihilation models, for the moment, please start with just the annihilation mode

Ratio

The flux of gamma-rays from DM in the case of annihilation is $$\phi(>E0)=\frac{1}{4\pi}\frac{\sigma v}{2m^2}\int{E_0}^m \frac{dN}{dE}dE\times J$$ where $m$ is the DM mass, $dN/dE$ is the photon spectrum per annihilation, $J$ is the astrophysical factor. Please tell me if you are not familiar with this equation.

From that equation, by resolving for $\sigma v$ one has

$$ \sigma v = \frac{8\pi m^2}{J} \frac{1}{\int_{E_0}^m \frac{dN}{dE}dE} $$

Now consider two cases: $b\bar{b}$ and $\tau^+\tau^-$. The only thing that change in the above equation for a given mass and target is the integral over the spectrum, from the energy threshold to the mass. This is the ratio that you have to compute, specifically:

$$ r = \frac{\int_{E0}^m \frac{dN (\tau^+\tau^-)}{dE}dE}{\int{E_0}^m \frac{dN (b\bar{b})}{dE}dE} $$

$r$ is computed solely from HDMS.

Application of casting

Once you have computed $r$ for these two channels and for the annihilation case you should

kumark0 commented 1 year ago

Thank you for clarifying, sir. I have a doubt regarding the 'r' that you mentioned. So for each dark matter mass, I have the 'r' values ranging from 0.1 to 10 (approx.) for all energy fractions x from the HDM code :

 dNdx_bb = HDMSpectra.spec(finalstate, initialstate1, x, m, data,annihilation=True)
dNdx_tautau = HDMSpectra.spec(finalstate, initialstate2, x, m, data,annihilation=True)
r=dNdx_bb/dNdx_tautau

To get the limit for each DM mass, I need to multiply only one ratio value, but HDM gives the ratio over the entire energy fraction x.

image

The above fig shows the range of limit for all the DM masses ranging from 503 GeV to 10TeV. As I had the 'r' for all x, once i multiplied it corresponding limit of bb , I got the range of values of sigmaV for each DM mass.

micheledoro commented 1 year ago

There is something wrong with the plot. The plot should have

The ratio

$$ r = \frac{\int_{E0}^m \frac{dN (\tau^+\tau^-)}{dE}dE}{\int{E_0}^m \frac{dN (b\bar{b})}{dE}dE} $$

should be computed from the integral of the spectrum not the spectrum itself. That's why your method .spec is not correct.

kumark0 commented 1 year ago

Hi @kumark0 apologies for the late reply. Comments

Annihilation/decay

There is no need to use both decay and annihilation models, for the moment, please start with just the annihilation mode

Ratio

The flux of gamma-rays from DM in the case of annihilation is ϕ(>E0)=14πσv2m2∫E0mdNdEdE×J where m is the DM mass, dN/dE is the photon spectrum per annihilation, J is the astrophysical factor. Please tell me if you are not familiar with this equation.

From that equation, by resolving for σv one has

σv=8πm2J1∫E0mdNdEdE

Now consider two cases: bb¯ and τ+τ−. The only thing that change in the above equation for a given mass and target is the integral over the spectrum, from the energy threshold to the mass. This is the ratio that you have to compute, specifically:

r=∫E0mdN(τ+τ−)dEdE∫E0mdN(bb¯)dEdE

r is computed solely from HDMS.

Application of casting

Once you have computed r for these two channels and for the annihilation case you should

  • take Segue bb¯ limit, multiply each mass point with r and compare with the Segue τ+τ− limit
  • take Segue τ+τ− limit, multiply each mass point with 1/r and compare with the Segue bb¯ limit

Sir, if we have defined the r= ration of tau by bb then shouldn't we need to multiply the bb llimmit with (1/r) to get the tautau limits and vice versa ?

kumark0 commented 1 year ago

Here is the output for the same. image

Code

`

cm = plt.get_cmap('seismic')

Specify the particles

finalstate = 22 # pdg id of the photon initialstate1 = 5 # pdg id of the b-quark initialstate2 = 15 #pdg id of tau leptops

dark matter mass in GeV, here 1 EeV

mDM =np.array([0.503618e3 ,0.727038e3 ,1.008575e3 ,1.325480e3 ,1.865110e3,2.762393e3 ,3.821204e3 ,5.346373e3 , 6.887659e3,10.0000e3])

x = np.logspace(-4.,0.,1000) # Energy fraction values, x = E/mDM

published limits on bb channel

bb_limits=np.array([8.249703e-24,6.542646e-24,5.812775e-24,5.544158e-24,5.388983e-24,5.703801e-24,6.094407e-24,6.957700e-24, 7.868472e-24,9.735452e-24]) tautau_limits=np.array([1.1654667955715847e-24, 1.1914159990017953e-24, 1.2180511137466625e-24, 1.2995040667614524e-24, 1.3862397753629218e-24, 1.4785458541736626e-24, 1.611204005705647e-24, 1.755504701987632e-24, 1.872232916109887e-24, 2.1281039743243496e-24 , 2.5784119955882113e-24, 2.993389123568225e-24, 3.402685818812826e-24, 3.786503459534035e-24, 4.396825205113306e-24, 4.998758394327381e-24, 5.806021997715327e-24, 6.3256394260508436e-24, 6.892168661563792e-24, 7.35153530261551e-24, 8.183683111871088e-24, 8.918731753954684e-24])

r=[] tautau_HDM_limit=[] bb_HDM_limit=[] plt.figure(figsize=(8,6))

Extract the spectrum using HDMSpectra.spec

for i, m in enumerate(mDM):

dNdx_bb = HDMSpectra.spec(finalstate, initialstate1, x, m, data,annihilation=True)
dNdE_bb=dNdx_bb/m
dNdx_tautau = HDMSpectra.spec(finalstate, initialstate2, x, m, data,annihilation=True)
dNdE_tautau=dNdx_tautau/m

Evals=x*m
r_bb = np.trapz(dNdE_bb, Evals)
r_tautau=np.trapz(dNdE_tautau,Evals)
ratio=r_tautau/r_bb
r.append(ratio)

for i in range(10):

bb_HDM_limit.append((1/r[i])*tautau_limits[i])

tautau_HDM_limit.append(bb_limits[i]*(1/r[i]))

`

micheledoro commented 1 year ago

Hi @kumark0, now that the comparison with HDMS is almost complete (see #22 ), do you wish to continue with this project of 'casting limits' or do you want to stop it?

kumark0 commented 1 year ago

Dear Prof. @micheledoro,

It would have been wonderful to continue the project, but since my master's thesis project has started, I won't be able to give the required time for this project. I'm looking forward to working with you in the future.

Also, I have been applying for a PhD. Graduate studies in the USA, and I would need to submit a reference letter for the same. Can you please send the reference letter? I would be glad to send the details required if you can write a letter of recommendation for my application.

With regards, Krishna