AtomDB / pyatomdb

AtomDB project (www.atomdb.org)
Docs: atomdb.readthedocs.io
Other
32 stars 12 forks source link

Building a cooling curve with AtomDB #16

Closed TarasKuzyo closed 5 years ago

TarasKuzyo commented 5 years ago

It is rather a feature request on the existing documentation.

I am trying to recreate a cooling curve for optically thin plasma (assuming CIE) $\Lambda(T)$ using pyatomdb. I was able to get line emissions by following the example in make_line_list.py, yet I am failing at getting the continuum component (I understand that it comes from apec_coco.fits but cannot find the right tool/method to process it).

Could you guide me in the right direction (or point to the appropriate resource)?

jagophile commented 5 years ago

Hi, I've added an example called "calc_power.py" to the examples directory, which should hopefully do what you need. Let me know if it's what you are looking for.

TarasKuzyo commented 5 years ago

@jagophile, Thank you -- it works great!

One last question: does ATOMDB data implies AG89 abundances by default, or do I need to re-scale separate line contributions when calculating the total cooling power?

jagophile commented 5 years ago

AG89 abundances are in there by default. So those powers are assuming for every 1 part hydrogen, there are 0.0977 parts helium , 0.000851 parts oxygen and so on (all of this is by number, not mass). If your plasma has a different composition you need to adjust accordingly.

Marking this as closed, though if you still have questions feel free to keep commenting

TarasKuzyo commented 5 years ago

@jagophile, I was wondering how do I get a cooling curve for NEI plasma for a given tau. I didn't find a corresponding pyatomdb.spectrum.make_spectrum function for NEI. Is there a way to adapt the CIE example tho this case?

jagophile commented 5 years ago

Yes - will post probably tomorrow

jagophile commented 5 years ago

Hello! This took a little longer than I hoped - sorry. I've added "calc_power_nei.py" to the examples directory.

TarasKuzyo commented 5 years ago

Thank you a lot! Did you push the example code to the master branch?

jagophile commented 5 years ago

Ah, I had put it in the wrong place. Should now be in the examples folder.

TarasKuzyo commented 5 years ago

@jagophile, I have a question regarding init_pop values.

A little background on the problem: I am carrying 3D simulations of a supernova remnant (Tycho SNR) and want to compare modeled thermal X-ray radiation with the observed map.

I am not sure what are the right values for init_pop array in pyatomdb.apec.solve_ionbal_eigen calculations. I came up with the following approach.

if Z > 2:
    init_pop[-1] = 0.5
    init_pop[-2] = 0.25
    init_pop[-3] = 0.25
else:
    init_pop[-1] = 1

Could you advice is those inital values are realistic or if there is a better way to set those up?

Also could you suggest a proper way to acknowledge your help with the cooling code? Is it ok to include you in the authors list or is there another way you prefer?

jagophile commented 5 years ago

Hi, Not sure what you are looking for here. init_pop is the population before ionizing starts. So typically, for an SNR, you would assume it is all neutral. So init_pop[0] = 1, everything else is 0.

The tau and temperature are used within the calculation, which calculates the real ionization population (takes the all neutral plasma, evolves it to something with n_e*t = 10^12, or whatever is appropriate for your analysis) and then calculated the appropriate spectra.

unless you have reason to believe that the plasma is strongly recombining, I wouldn't start by populating the fully stripped states.

But I don't know what your model is, so I can't immediately help without more information. I'd suggest email (afoster@cfa.harvard.edu) with more details.