Non-parallel rays for source broadening

[lemmatum]

I've been making a simulation for a Johann spectrometer based on the example in the repository, and I noticed that the rays from the GeometricSource to the crystal are parallel. Like so:

When I reduce this to a simple example of a GeometricSource and a detector it confirms that the rays do not diverge, and the size of the source is simply imprinted directly onto the detector

Is there a way to simulate diverging (non-parallel) rays in order to simulate line broadening due to the finite size of the source?

[kklmn]

I've run the example examples/withRaycing/06_AnalyzerBent1D/01B_SourceZCrystalThetaAlpha.py out of the box (except making showIn3D = True and isJohansson = False # i.e. is Johann and do see a divergent beam:

A similar picture I see for the example of 2D bent crystal examples/withRaycing/07_AnalyzerBent2D/01BD_SourceZCrystalThetaAlpha.py. So you're probably using an example for something else or have modified the source yourself.

Linear source distributions of GeometricSource (RTFM) are defined by distx, disty and distz parameters, whereas angular distributions by distxprime and distzprime. The distribution sizes are correspondingly given by dx, dy, dz, dxprime and dzprime.

[lemmatum]

Yes, I updated the geometry to my specific problem, and I forgot to update the value for dxprime and dzprime. When I increase the value to 0.07 then I do see a diverging source (thanks!):

What is the physical meaning of dxprime and dzprime? Would this be something like the solid angle subtended by the Bragg crystal relative to the source?

From examples/withRaycing/06_AnalyzerBent1D/01B_SourceZCrystalThetaAlpha.py I see:

pdp = 2. * R * math.sin(theta + alpha - dyCrystal/6/R)
beamLine.source.dxprime = dxCrystal / pdp
beamLine.source.dzprime = dyCrystal * math.sin(theta - alpha) / pdp

So it appears to be something like the size of the crystal divided by the distance traversed between source and detector on the Rowland Circle. This looks almost like a solid angle, but not quite.

[lemmatum]

As a follow-up question. Is there a way to extract the plotted data into numpy arrays for further analysis?

If I defined plotDetE = xrtp.XYCPlot(...) and then try to access the underlying data via plotDetE.yaxis.data it asks for a beam, and I'm unsure what the method wants.

I am mainly interested in this for simulating the dispersion of the spectrometer and then fitting the dispersion curve with some function.

[kklmn]

What is the physical meaning of dxprime and dzprime?

Please tell what is unclear in the above-cited docs.

Is there a way to extract the plotted data into numpy arrays for further analysis?

Saving the results

[lemmatum]

Please tell what is unclear in the above-cited docs.

For example, the docs for GeometricSource just says

dx, dy, dz, dxprime, dzprime: float
Linear (dx, dy, dz) and angular (dxprime, dzprime) source sizes. For normal distribution is sigma or (sigma, cut_limit), for flat is full width or tuple (min, max), for annulus is tuple (rMin, rMax), otherwise is ignored.

So from the docs, all I know is that dxprime and dzprime somehow relate to angular source size, but it doesn't tell me anything about the units (mm, radians, degrees), or how I would make sure to fill an object with rays. After playing around with my own minimal example of a point source and a detector (no crystal) it looks like dxprime is the angle in radians subtended by the detector out of 2 * pi. It would help to clarify this in the docstring for GeometricSource.

Saving the results

Thanks! I'll take a look at implementing an afterScript.

[kklmn]

Units

or how I would make sure to fill an object with rays

The workflow is outlined right in the next section under Units.

[lemmatum]

Okay, so I can now extract things like intensity and angle of rays incident at the detector. Is there a way to extract the photon energy (or energies) of the rays which hit a specific "pixel" at the detector plane? Would this be what is called the 4D mutual intensity?

[kklmn]

To analyze beams on a specific plane, we create a screen there and use its pixels for histogramming the beams. All rays coming to a particular pixel lose their identity and are averaged over the pixel. Energy is averaged as well in the sense of color mixture. If you don't want histogramming, you can do your own analysis at the end of run_process().

Would this be what is called the 4D mutual intensity?

What this? Please ask again.

[lemmatum]

Extracting the averaged photon energy as binned from the histogram should be fine for my purposes. Is there a specific method or attribute I can use in run_process to extract these photon energies as either a 1D or 2D spatial array?

Oh, I was wondering whether the method for 4D mutual intensity is a method for extracting either the individual photon energies of the rays or the average (histogram binned) energy or whether it does something else entirely.

[kklmn]

run_process() is designed to start in parallel in several Python contexts. This method doesn't do histogramming, it feeds histogramming. Any beam container has an attribute .E -- an energy array per ray, also .x, .y and .z among others. You can analyze them the way you want right in run_process().
After histogramming, you can extract your energy from the plot container where the histogramming was done. For doing this, you need to define the caxis as energy, which is a default plot setting.

Mutual intensity is used in coherence applications, this is not your case.

[lemmatum]

plot.yaxis.total1D gives the flux per pixel in the y-direction, but instead I want the photon energy.
plot.caxis.binCenters gives the photon energies for the caxis bins, but I don't see how to connect these photon energies to what is on the yaxis. Do any of the plot.yaxis methods or attributes have information on photon energy, or some index to connect it to the binning on the caxis?

There is maybe some connection through plot.yaxis.total1D_RGB and plot.caxis.total1D_RGB but this seems like a very round about way of trying to obtain photon energy values for the yaxis.

[lemmatum]

I got something close to what I want by using raycing.get_energy and raycing.get_z under define_plots:

Would there be a way for me to generate my own matplotlib plot under define_plots? The above data would be easier to visualize as a regular scatter plot instead of an XYCPlot. Alternatively, I could extract the underlying data of the energy and z axes in the afterScript.

Do you have any suggestions on what is the most simple way to do this using your code?

[kklmn]

You're trying to solve an undefined problem. The color mixing in our plots is done in the sense of RGB addition and is good for visual representation only. In general, you have a whole spectrum per pixel. total1D_RGB will give you intensity weighted averaged energy but the width of energy distribution will be lost. If you see a cyan color, you don't know how it was obtained: a single color, a sum of blue and green, a sum of red and violet, and any combination of those.

I suggest studying monochromatic cases at a scanned single photon energy.

[lemmatum]

I agree, the RGB method would be too poorly defined.

Is there a similar mixing for raycing.get_energy or can I consider each element in that array as monochromatic instead of a blend of energies? For the plot in my previous comment it seems to give a distinct energy for each given line and for each pixel, which is what I was expecting (and wanting).