import SmartPhantom as SmrtPhntm

import numpy as np

"""  --------    --------    --------    -------- 

  Script:

     --------    --------    --------    -------- """

#--------  --------  Initialisation  --------  -------- 
# Set debug flag:
Debug = False
print("Starting: debug = ", Debug)

# Set plot flag:
Plot = True
print("          plot = ", Plot)

# Define water box:
Edges = SmrtPhntm.defineWaterBox(300., 300., 300.)
Vox   = SmrtPhntm.defineVoxel(1., 1., 1.)
print(" Water box width (+/-):", Edges[0], " mm, height (+/-):", Edges[1], " mm, length (+/-):", Edges[2], " mm")
print(" Voxel size (x, y, z):", Vox, " mm")

# Define time of flight max and time step:
ToFMax = 100.
TStep  = 1.
print(" Max time of flight and time step:", ToFMax, TStep)

# Define array to store energy sum:
NBinsX, NBinsY, NBinsZ = SmrtPhntm.getNumberOfVoxels(Edges, Vox)
if Debug:
     print(" NBinsX, NBinsY, NBinsZ:", NBinsX, NBinsY, NBinsZ)
NBinsT = SmrtPhntm.getNumberOfTimeSteps(ToFMax, TStep)
if Debug:
     print(" NBinsT:", NBinsT)
     
#--------  --------  Initialisation -- end  --------  -------- 

#--------  --------  Load power spectrum --------  -------- 
    
Power = np.load('12-PowerSpectra/power-spectrum-protons_only.npy')

print(" Power array loaded.")

#--------  --------  Load power spectrum -- end  --------  -------- 
 
#--------  --------  Make plots, to check(!) --------  -------- 
    
if Plot:
     Dummy = SmrtPhntm.plotPower3Coord(NBinsX, NBinsY, NBinsZ, Edges, Power)

#--------  --------  Make plots, to check(!) -- end  --------  -------- 

print("Complete")
