diff --git a/ess/lib_tw.py b/ess/lib_tw.py
index 33279c32c7c772bc4d26f972352a384b3cd581d4..43c32672680bb7de682248565ecc46e934c006db 100644
--- a/ess/lib_tw.py
+++ b/ess/lib_tw.py
@@ -429,7 +429,6 @@ class PROJECT:
 
 # ---- Data related
 
-
 class PARTRAN:
     """
         Note:
@@ -438,143 +437,155 @@ class PARTRAN:
 
         History:
 
-        - 2016.02.17: Changed how to identify the line of indices.
-        - 2016.02.17: Added a logic to avoid #/0 for LEBT.
+        - 2016-02-17: Changed how to identify the line of indices.
+        - 2016-02-17: Added a logic to avoid #/0 for LEBT.
+        - 2020-08-28: Adopted to python3.
+        - 2020-08-28: 99% emittance added.
+        - 2020-08-28: alfz was wrong (in power of gamma) and fixed.
 
     """
 
     def __init__(self, file_name):
+        
+        # Consts for conversions
+        mc2  = 938.272088              # MeV
+        wlen = 2.99792e8/352.21e6*1e3  # mm
 
-        # Consts to convert phs to z.
-        c = 2.99792458
-        freq = 352.21
-
-        # Extract indices.
+        # Extract indices
         with open(file_name) as file:
+            i_header = 0
             for lin in file.readlines():
                 lin = lin.split()
                 if "##" in lin[0]:
-                    idx_s = lin.index("z(m)")
+                    idx_Nelem = lin.index("##")
+                    idx_s     = lin.index("z(m)")
                     idx_gamma = lin.index("gama-1")
-                    idx_x = lin.index("x0")
-                    idx_y = lin.index("y0")
-                    idx_phs = lin.index("p0")
-                    idx_sigx = lin.index("SizeX")
-                    idx_sigy = lin.index("SizeY")
-                    idx_sigz = lin.index("SizeZ")
-                    idx_sigp = lin.index("SizeP")
-                    idx_alfx = lin.index("sxx'")
-                    idx_alfy = lin.index("syy'")
-                    idx_alfz = lin.index("szdp")
-                    idx_epsx = lin.index("ex")
-                    idx_epsy = lin.index("ey")
-                    idx_epsz = lin.index("ezdp")
-                    idx_epsp = lin.index("ep")
-                    idx_halx = lin.index("hx")
-                    idx_haly = lin.index("hy")
-                    idx_halp = lin.index("hp")
+                    idx_x     = lin.index("x0")
+                    idx_y     = lin.index("y0")
+                    idx_phs   = lin.index("p0")
+                    idx_sigx  = lin.index("SizeX")
+                    idx_sigy  = lin.index("SizeY")
+                    idx_sigz  = lin.index("SizeZ")
+                    idx_sigp  = lin.index("SizeP")
+                    idx_alfx  = lin.index("sxx'")
+                    idx_alfy  = lin.index("syy'")
+                    idx_alfz  = lin.index("szdp")  # Note unit is [um]
+                    idx_epsx  = lin.index("ex")
+                    idx_epsy  = lin.index("ey")
+                    idx_epsz  = lin.index("ezdp")
+                    idx_epsp  = lin.index("ep")
+                    idx_e99x  = lin.index("ex99")
+                    idx_e99y  = lin.index("ey99")
+                    idx_e99z  = lin.index("ep99")  # Note 'ep99' is actually for z-z'
+                    idx_halx  = lin.index("hx")
+                    idx_haly  = lin.index("hy")
+                    idx_halp  = lin.index("hp")
                     idx_Nptcl = lin.index("npart")
-                    idx_loss = lin.index("Powlost")
+                    idx_loss  = lin.index("Powlost")
                     break
-
-        # Extract data.
-        with open(file_name) as file:
-            data = []
-            flag = 0
-            for lin in file.readlines():
-                lin = lin.split()
-                if flag == 1:
-                    data.append(map(float, lin))
-                if "##" in lin[0]:
-                    flag = 1
-            data = numpy.array(data).transpose()
-
+                i_header += 1
+                                        
+        # Extract data
+        data = numpy.loadtxt(file_name, skiprows=i_header+1).transpose()
+            
         # Instances
-        self.s = data[idx_s]
-        self.x = data[idx_x]
-        self.y = data[idx_y]
-        self.phs = data[idx_phs]
-        self.sigx = data[idx_sigx]
-        self.sigy = data[idx_sigy]
-        self.sigz = data[idx_sigz]
-        self.sigp = data[idx_sigp]
-        self.epsx = data[idx_epsx]
-        self.epsy = data[idx_epsy]
-        self.epsz = data[idx_epsz]
-        self.epsp = data[idx_epsp]
-        self.halx = data[idx_halx]
-        self.haly = data[idx_haly]
-        self.halz = data[idx_halp]
-        self.halp = data[idx_halp]
+        self.Nelem = data[idx_Nelem].astype(numpy.int)
+        self.s     = data[idx_s]
+        self.x     = data[idx_x]
+        self.y     = data[idx_y]
+        self.phs   = data[idx_phs]
+        self.sigx  = data[idx_sigx]
+        self.sigy  = data[idx_sigy]
+        self.sigz  = data[idx_sigz]
+        self.sigp  = data[idx_sigp]
+        self.epsx  = data[idx_epsx]
+        self.epsy  = data[idx_epsy]
+        self.epsz  = data[idx_epsz]
+        self.epsp  = data[idx_epsp]
+        self.e99x  = data[idx_e99x]
+        self.e99y  = data[idx_e99y]
+        self.e99z  = data[idx_e99z]
+        self.halx  = data[idx_halx]
+        self.haly  = data[idx_haly]
+        self.halz  = data[idx_halp]
+        self.halp  = data[idx_halp]
         self.Nptcl = data[idx_Nptcl]
-        self.loss = data[idx_loss]
+        self.loss  = data[idx_loss]
 
         # Way around for 0 emitt
-        for i in range(len(self.s)):
-            if self.epsx[i] == 0.0:
-                self.epsx[i] = numpy.inf
-            if self.epsy[i] == 0.0:
-                self.epsy[i] = numpy.inf
-            if self.epsz[i] == 0.0:
-                self.epsz[i] = numpy.inf
-            if self.epsp[i] == 0.0:
-                self.epsp[i] = numpy.inf
+        for i in numpy.arange(len(self.s)):
+            if self.epsx[i] == 0.0: self.epsx[i] = numpy.inf
+            if self.epsy[i] == 0.0: self.epsy[i] = numpy.inf
+            if self.epsz[i] == 0.0: self.epsz[i] = numpy.inf
+            if self.epsp[i] == 0.0: self.epsp[i] = numpy.inf
+            if self.e99x[i] == 0.0: self.e99x[i] = numpy.inf
+            if self.e99y[i] == 0.0: self.e99y[i] = numpy.inf
+            if self.e99z[i] == 0.0: self.e99z[i] = numpy.inf
 
         # Additional instances
         self.gamma = data[idx_gamma] + 1.0
-        self.beta = numpy.sqrt(1.0 - 1.0 / self.gamma ** 2)
-        self.z = -self.phs * self.beta * (c / freq * 1e5) / 360.0
-        self.betx = self.sigx ** 2 / self.epsx * self.beta * self.gamma
-        self.bety = self.sigy ** 2 / self.epsy * self.beta * self.gamma
-        self.betz = self.sigz ** 2 / self.epsz * self.beta * self.gamma ** 3
-        self.betp = self.sigp ** 2 / self.epsp
-        self.alfx = -data[idx_alfx] / self.epsx * self.beta * self.gamma
-        self.alfy = -data[idx_alfy] / self.epsy * self.beta * self.gamma
-        self.alfz = -data[idx_alfz] / self.epsz * self.beta * self.gamma ** 3
-        self.alfp = -self.alfz
-
+        self.beta  = numpy.sqrt(1.0 - 1.0 / self.gamma ** 2)
+        self.z     = -self.phs * self.beta * wlen / 360.0
+        self.betx  = self.sigx ** 2 / self.epsx * self.beta * self.gamma
+        self.bety  = self.sigy ** 2 / self.epsy * self.beta * self.gamma
+        self.betz  = self.sigz ** 2 / self.epsz * self.beta * self.gamma ** 3
+        self.betp  = self.sigp ** 2 / self.epsp
+        self.alfx  = -data[idx_alfx] / self.epsx * self.beta * self.gamma
+        self.alfy  = -data[idx_alfy] / self.epsy * self.beta * self.gamma
+        self.alfz  = -data[idx_alfz] / self.epsz * self.beta * self.gamma
+        self.alfp  = -self.alfz
+        self.e99p  = self.e99z * 1e-3 * 360.0 * mc2 / wlen
+        
+        print(data[idx_alfz][ 0])
+        print(data[idx_alfz][-1])
+        
         # Set inf emitt back to 0
-        for i in range(len(self.s)):
-            if self.epsx[i] == numpy.inf:
-                self.epsx[i] = 0.0
-            if self.epsy[i] == numpy.inf:
-                self.epsy[i] = 0.0
-            if self.epsz[i] == numpy.inf:
-                self.epsz[i] = 0.0
-            if self.epsp[i] == numpy.inf:
-                self.epsp[i] = 0.0
+        for i in numpy.arange(len(self.s)):
+            if self.epsx[i] == numpy.inf: self.epsx[i] = 0.0
+            if self.epsy[i] == numpy.inf: self.epsy[i] = 0.0
+            if self.epsz[i] == numpy.inf: self.epsz[i] = 0.0
+            if self.epsp[i] == numpy.inf: self.epsp[i] = 0.0
+            if self.e99x[i] == numpy.inf: self.e99x[i] = 0.0
+            if self.e99y[i] == numpy.inf: self.e99y[i] = 0.0
+            if self.e99z[i] == numpy.inf: self.e99z[i] = 0.0
+            if self.e99p[i] == numpy.inf: self.e99p[i] = 0.0
 
         # Convert to list (not necessary?)
-        self.s = self.s.tolist()
+        self.Nelem = self.Nelem.tolist()
+        self.s     = self.s.tolist()
         self.gamma = self.gamma.tolist()
-        self.beta = self.beta.tolist()
-        self.x = self.x.tolist()
-        self.y = self.y.tolist()
-        self.z = self.z.tolist()
-        self.phs = self.phs.tolist()
-        self.sigx = self.sigx.tolist()
-        self.sigy = self.sigy.tolist()
-        self.sigz = self.sigz.tolist()
-        self.sigp = self.sigp.tolist()
-        self.betx = self.betx.tolist()
-        self.bety = self.bety.tolist()
-        self.betz = self.betz.tolist()
-        self.betp = self.betp.tolist()
-        self.alfx = self.alfx.tolist()
-        self.alfy = self.alfy.tolist()
-        self.alfz = self.alfz.tolist()
-        self.alfp = self.alfp.tolist()
-        self.epsx = self.epsx.tolist()
-        self.epsy = self.epsy.tolist()
-        self.epsz = self.epsz.tolist()
-        self.epsp = self.epsp.tolist()
-        self.halx = self.halx.tolist()
-        self.haly = self.haly.tolist()
-        self.halz = self.halz.tolist()
-        self.halp = self.halp.tolist()
+        self.beta  = self.beta.tolist()
+        self.x     = self.x.tolist()
+        self.y     = self.y.tolist()
+        self.z     = self.z.tolist()
+        self.phs   = self.phs.tolist()
+        self.sigx  = self.sigx.tolist()
+        self.sigy  = self.sigy.tolist()
+        self.sigz  = self.sigz.tolist()
+        self.sigp  = self.sigp.tolist()
+        self.betx  = self.betx.tolist()
+        self.bety  = self.bety.tolist()
+        self.betz  = self.betz.tolist()
+        self.betp  = self.betp.tolist()
+        self.alfx  = self.alfx.tolist()
+        self.alfy  = self.alfy.tolist()
+        self.alfz  = self.alfz.tolist()
+        self.alfp  = self.alfp.tolist()
+        self.epsx  = self.epsx.tolist()
+        self.epsy  = self.epsy.tolist()
+        self.epsz  = self.epsz.tolist()
+        self.epsp  = self.epsp.tolist()
+        self.e99x  = self.e99x.tolist()
+        self.e99y  = self.e99y.tolist()
+        self.e99z  = self.e99z.tolist()
+        self.e99p  = self.e99p.tolist()
+        self.halx  = self.halx.tolist()
+        self.haly  = self.haly.tolist()
+        self.halz  = self.halz.tolist()
+        self.halp  = self.halp.tolist()
         self.Nptcl = self.Nptcl.tolist()
-        self.loss = self.loss.tolist()
-
+        self.loss  = self.loss.tolist()
+        
     def loss_den(self, file_name_dt="", dlt_dt=5e-6):
 
         return loss_elem2den(self.s, self.loss, file_name_dt, dlt_dt)