content: /*
   * Copyright (C) 2014 - Andreas Maier 
   * CONRAD is developed as an Open Source project under the GNU General Public License (GPL).
   */
  package edu.stanford.rsl.tutorial.physics;
  
  import ij.gui.Plot;
  
  import java.util.ArrayList;
  
  import edu.stanford.rsl.conrad.geometry.shapes.simple.Edge;
  import edu.stanford.rsl.conrad.geometry.shapes.simple.PointND;
  import edu.stanford.rsl.conrad.geometry.trajectories.Trajectory;
  import edu.stanford.rsl.conrad.physics.PhysicalObject;
  import edu.stanford.rsl.conrad.physics.PolychromaticXRaySpectrum;
  import edu.stanford.rsl.conrad.physics.absorption.PolychromaticAbsorptionModel;
  import edu.stanford.rsl.conrad.physics.absorption.SelectableEnergyMonochromaticAbsorptionModel;
  import edu.stanford.rsl.conrad.physics.materials.database.MaterialsDB;
  import edu.stanford.rsl.conrad.utils.Configuration;
  
  /**
   * Example class to play and test spectral absorption. Here we generate a spectrum and plot it.
   * Next, we create a monochromatic and a polychromatic absorption model. We investigate both in terms of their absorption behavior and create a plot for the mono- and the polychromatic case.
   * The difference can be looked at, both both are plotted in the same coordinate system. This code is also featured on the CONRAD tutorial website in more detail.
   * @author Andreas Maier
   *
   */
  public class SpectralAbsorption {
  
  	public static void main(String [] args){
  
  		// Enable access to the global configuration
  		Configuration.loadConfiguration();
  
  		// Fetch some geometry information from the global configuration
  		Trajectory trajectory = Configuration.getGlobalConfiguration().getGeometry();
  		double detectorArea = trajectory.getPixelDimensionX() * trajectory.getDetectorWidth() * trajectory.getDetectorHeight() * trajectory.getPixelDimensionY();
  		double pixelArea = trajectory.getPixelDimensionX() * trajectory.getPixelDimensionY();
  
  		// set spectral parameters
  		double min = 10;
  		double max = 60;
  		double delta = 0.5;
  		double peakVoltage = 40;			
  		double timeCurrentProduct = 1;	
  
  		// create spectrum
  		PolychromaticXRaySpectrum spectrum = new PolychromaticXRaySpectrum(min, max, delta, peakVoltage, timeCurrentProduct);
  		// get total photon flux in [photons/mm^-2]
  		double totalPhotonFlux = spectrum.getTotalPhotonFlux();
  		/* input energy in J */
  		double inputEnergy = peakVoltage * timeCurrentProduct;
  		/* convert to keV */
  		inputEnergy *= 6.24150934E15; 
  		/* output energy in keV */
  		double outputEnergy = detectorArea * totalPhotonFlux * spectrum.getAveragePhotonEnergy();
  		/* compute source efficiency, i.e. the amount of energy that is detected compared to the amount of energy that is used for creation of the spectrum */
  		double ratio = outputEnergy / inputEnergy;
  		// create spectral plot
  		Plot plot = new Plot("C-Arm Spectrum with "  + detectorArea * totalPhotonFlux + " total photons and " + pixelArea * totalPhotonFlux + " photons per pixel, avg energy "+ spectrum.getAveragePhotonEnergy() + " [keV] and " +ratio*100+ " % energy efficiency", "Energy [keV]", "Photon Flux [photons/mm^2]", spectrum.getPhotonEnergies(), spectrum.getPhotonFlux());
  		plot.show();
  
  		// create absorption / length plots
  		// initialize arrays
  		int vals = 20;
  		double lengths [] = new double [vals];
  		double monoAttenuation []= new double[vals];
  		double polyAttenuation []= new double[vals];
  		// create monochromatic absorption model
  		SelectableEnergyMonochromaticAbsorptionModel monochromaticAbsorptionModel = new SelectableEnergyMonochromaticAbsorptionModel();
  		monochromaticAbsorptionModel.configure(spectrum.getAveragePhotonEnergy());
  		// create polychromatic absorption model
  		PolychromaticAbsorptionModel polychromaticAbsorptionModel = new PolychromaticAbsorptionModel();
  		polychromaticAbsorptionModel.setInputSpectrum(spectrum);
  		try {
  			polychromaticAbsorptionModel.configure();
  		} catch (Exception e) {
  			// TODO Auto-generated catch block
  			e.printStackTrace();
  		}
  		// create list of line intersections.
  		// Here we have a line with only one intersection
  		// of water
  		PhysicalObject segment = new PhysicalObject();
  		segment.setMaterial(MaterialsDB.getMaterial("water"));
  		segment.setNameString("Water Path");
  		ArrayList<PhysicalObject> segments = new ArrayList<PhysicalObject>();
  		segments.add(segment);
  		// Iterate over the different intersection lengths
  		for (int i = 0; i < vals; i++){
  			lengths[i] = i*20;
  			// create intersection path
  			Edge edge = new Edge(new PointND(0), new PointND(lengths[i]));
  			segment.setShape(edge);
  			// evaluate attenuation models
  			monoAttenuation[i]= monochromaticAbsorptionModel.evaluateLineIntegral(segments);
  			polyAttenuation[i]= polychromaticAbsorptionModel.evaluateLineIntegral(segments);
  		}
  		// create plots
  		Plot plot2 = new Plot("Monochromatic absorption q, avg energy "+ spectrum.getAveragePhotonEnergy() + " [keV]", "path length in water", "q_mono", lengths, monoAttenuation);
  		plot2.show();
  		Plot plot3 = new Plot("Polychromatic absorption q, avg energy "+ spectrum.getAveragePhotonEnergy() + " [keV]", "path length in water", "q_poly", lengths, polyAttenuation);
  		plot3.show();
  	}
  
  
  }
  