Core code analysis: How GISAS intensity is computed

last updated 05dec25

ISimulation::simulate()
  prepareSimulation() virtual
  re_sample = ReSample::make(..)
  runSingleSimulation(re_sample, batch...)
    threaded...
      loop over i_pixel in thread
        runComputation(re_sample, i_pixel, wgt) virtual = 0

ScatteringSimulation::runComputation(re_sample, i_pixel, wgt)
  // compute intensity from scattering and specular reflection
  // multiply with beam intensity, footprint factor, and weight

Evaluate::scatteredContribution(re_sample, ele, options, polarized) -> double
  // implemented in Sim/Evaluate/EvaluateScattering.cpp
  sum over s: re_sample.structs()
    intensity += scatteringWeight(s) * scatteringIntensity(ele, options, lambda s->scalarScattering(ele))
    // or polarized: s->polarizedScattering(ele)
  // where scatteringWeight(s) dispatches to:
  //   s2d->areaNumberDensity() * s->coverage()   for Struct2D
  //   s3d->volumeNumberDensity() * s->coverage() for Struct3D

// For all 2D structures (Sample/Aggregate/Struct2D.cpp):
Struct2D::scalarScattering(ele) -> double
  diffuse_intensity, amplitude := 0
  loop(csp: coheringPiles())
     ff = csp->summedFF(ele)
     fraction = csp->relativeAbundance()
     amplitude += fraction * ff
     diffuse_intensity += fraction * |ff|^2
  q := ele.meanQ()
  return diffuse_intensity + (structureFactor(q) - 1) * |amplitude|^2
  // structureFactor() returns 1.0 in base class Struct2D (random assemblies)
  // structureFactor() is overridden in Ordered2D subclasses

// CoheringPile stores particles with fixed relative positions (Sample/Scattering/CoheringPile.cpp)
CoheringPile::summedFF(ele) -> complex
  sum over spr: m_particles_in_pile
    += spr->dwbaFF(ele)

// DWBA computation (Sample/Scattering/IScatterer.cpp)
// IScatterer is base class of IParticle and Struct3D
IScatterer::dwbaFF(ele) -> complex
  // For pure Born approximation (only_outgoing):
  return theFF(wavevectors)
  // For full DWBA, compute 4 terms (T_in, R_in) x (T_out, R_out):
  return T_in * T_out * theFF(q_TT)   // direct scattering
       + R_in * T_out * theFF(q_RT)   // reflection then scattering
       + T_in * R_out * theFF(q_TR)   // scattering then reflection
       + R_in * R_out * theFF(q_RR)   // reflection, scattering, reflection

IParticle::theFF(WavevectorInfo wavevectors) -> complex
  // applies position phase and rotation
  return e^(i*q*position) * formfactor->theFF(rotated wavevectors)

IFormfactor::theFF(WavevectorInfo wavevectors) -> complex
  returns formfactor(wavevectors.getQ())

IFormfactor::formfactor(C3 q) -> complex, virtual = 0

Structure factor

Paracrystal2D:

Ordered2D::structureFactor(q, outer_iff)
  returns iff_no_inner(q, outer_iff)

Ordered2D::iff_no_inner(q, outer_iff)
  returns DWfactor(q) * (iff_without_dw(q) * outer_iff - 1.0) + 1.0

ICrystal2D::iff_without_dw(q)
  performs integration or simply
  returns interferenceForXi(m_lattice->rotationAngle(), q.x(), q.y())

Paracrystal2D::interferenceForXi(xi, qx, qy)
  returns interference1D(qx, qy, xi, 0)
    * interference1D(qx, qy, xi + m_lattice->latticeAngle(), 1)

RadialParacrystal:

RadialParacrystal::iff_without_dw(q)
  directly performs Hosemann computation