import sys, os
sys.path.insert(0, "/Users/kaestner/git/lib/lib")
sys.path.insert(0, "/Users/kaestner/git/scripts/python/")
if 'LD_LIBRARY_PATH' not in os.environ:
    os.environ['LD_LIBRARY_PATH'] = '/Users/kaestner/git/lib/lib'
   # os.execv(sys.argv[0], sys.argv)


import numpy as np
import muhrectomo as mt
import matplotlib.pyplot as plt
import amglib.readers as io
import amglib.imageutils as amg

import importlib
importlib.reload(mt)
importlib.reload(io)
importlib.reload(amg)

"""
Explanation: MuhRec in python demonstration
End of explanation
"""


recon = mt.Reconstructor(mt.bpMultiProj)

print("Created reconstructor :",recon.name())

"""
Explanation: Create a reconstructor object
End of explanation
"""


Nproj = 10
proj = np.ones([Nproj,256,256])

args = {"angles" : np.linspace(0,180,num=Nproj), 
        "weights" : np.ones(Nproj)}

recon.configure({   "center" : 50, 
                    "resolution" : 0.135
                })

recon.process(proj,args)

vol = recon.volume()

plt.imshow(vol[0])

"""
Explanation: Reconstruction workflow
Numerical dummy
This dummy is only to test the connection and that the reconstructor starts at all. Only nonsense data.
End of explanation
"""


ob = io.readImages('/data/P08062_wood/raw_CCD/ob_{0:04d}.fits',1,5,averageStack=True).mean(axis=0)
dc = io.readImages('/data/P08062_wood/raw_CCD/dark_{0:04d}.fits',1,5,averageStack=True).mean(axis=0)

proj = io.readImages('/data/P08062_wood/raw_CCD/wood_{0:04d}.fits',1,626) # This takes a while

fig,ax = plt.subplots(1,3,figsize=[15,7])
ax[0].imshow(dc)
ax[1].imshow(ob)
ax[2].imshow(proj[0])


"""
Explanation: The wood data
The wood data is the data I use a lot for demos, tests, and tutorials. It a piece of petrified wood DOI. 
Load projection data
End of explanation
"""


nproj = amg.normalizeImage(img=proj, ob=ob, dc=dc, doseROI=[100,250,200,300])

cproj = nproj[:,50:920,250:850]
plt.imshow(cproj[0])
del nproj

"""
Explanation: Preprocessing
Here the projections needs to be prepared for reconstruction. Typical operations are 
- Cropping
- Normalization (possibly with scattering correction)
- Spot cleaning
- Ring cleaning
End of explanation
"""


Nproj = cproj.shape[0]
# Information per projection
args = {"angles"  : np.linspace(0,360,num=Nproj), 
        "weights" : np.ones(Nproj)/Nproj}

# Geometry information
recon.configure({   "center" : 295, 
                    "resolution" : 0.05
                })

recon.process(cproj[:,500:600,:],args) # Reconstruct a part of the slices (32 slices here)

vol = recon.volume() # Retrieve the reconstructed volume

plt.imshow(vol[0])

import ipyvolume as ipv

ipv.quickvolshow(vol)

vol.shape

"""
Explanation: Prepare and run the back-projection
End of explanation
"""


