Skip to content
Snippets Groups Projects
TraceWin.py 46.1 KiB
Newer Older
Yngve Levinsen's avatar
Yngve Levinsen committed
        if self.version >= 6 and partran:
            coordinates += ["min_pos_moy", "max_pos_moy"]
            coordinate_units += ["m", "m"]
Yngve Levinsen's avatar
Yngve Levinsen committed
        for val, unit in zip(arrays, array_units):
            data_set = group.create_dataset(val, (length,), dtype="f")
Yngve Levinsen's avatar
Yngve Levinsen committed
            data_set[...] = getattr(self, val)
                data_set.attrs["unit"] = unit
Yngve Levinsen's avatar
Yngve Levinsen committed
        for val, unit in zip(coordinates, coordinate_units):
            data_set = group.create_dataset(val, (length, 7), dtype="f")
Yngve Levinsen's avatar
Yngve Levinsen committed
            data_set[...] = getattr(self, val)
                data_set.attrs["unit"] = unit
Yngve Levinsen's avatar
Yngve Levinsen committed
        if self.version >= 7 and partran:
            # 3 numbers per location..
            emit_data = ["rms_emit", "rms_emit2"]
            emit_units = ["m*rad", "m*m*rad*rad"]
Yngve Levinsen's avatar
Yngve Levinsen committed
            for val, unit in zip(emit_data, emit_units):
                data_set = group.create_dataset(val, (length, 3), dtype="f")
Yngve Levinsen's avatar
Yngve Levinsen committed
                data_set[...] = getattr(self, val)
                    data_set.attrs["unit"] = unit
        if partran:
Yngve Levinsen's avatar
Yngve Levinsen committed
            # 1 numbers per location and per run..
            data = ["lost", "powlost"]
            units = ["", "W"]
Yngve Levinsen's avatar
Yngve Levinsen committed
            for val, unit in zip(data, units):
                data_set = group.create_dataset(val, (length, self.Nrun), dtype="f")
Yngve Levinsen's avatar
Yngve Levinsen committed
                data_set[...] = getattr(self, val)
                    data_set.attrs["unit"] = unit
Yngve Levinsen's avatar
Yngve Levinsen committed
        fout.close()
class remote_data_merger:
    def __init__(self, base="."):
Yngve Levinsen's avatar
Yngve Levinsen committed
        self._base = base
        self._files = []
Yngve Levinsen's avatar
Yngve Levinsen committed
    def add_file(self, filepath):
        import os

        if os.path.exists(filepath):
Yngve Levinsen's avatar
Yngve Levinsen committed
            fname = filepath
Yngve Levinsen's avatar
Yngve Levinsen committed
            fullpath = os.path.join(self._base, filepath)
            if os.path.exists(fullpath):
Yngve Levinsen's avatar
Yngve Levinsen committed
                fname = fullpath
Yngve Levinsen's avatar
Yngve Levinsen committed
                raise ValueError("Could not find file " + filepath)
        if fname not in self._files:
            self._files.append(fname)

Yngve Levinsen's avatar
Yngve Levinsen committed
    def generate_partran_out(self, filename=None):
        Creates a string to be written to file
        each line is a list.

        If filename is given, writes directly to output file.

Yngve Levinsen's avatar
Yngve Levinsen committed

        import numpy as np

Yngve Levinsen's avatar
Yngve Levinsen committed
        h1 = []
        h2 = []
Yngve Levinsen's avatar
Yngve Levinsen committed

Yngve Levinsen's avatar
Yngve Levinsen committed
        d1 = []
        d2 = []
        d3 = []
Yngve Levinsen's avatar
Yngve Levinsen committed

        if self._files:
            for f in self._files:
                string = open(f, "r").read()
                split = string.split("$$$")
                if split[9] != "Data_Error":
                    raise ValueError("Magic problem, please complain to Yngve")

                thisdata = split[10].strip().split("\n")
Yngve Levinsen's avatar
Yngve Levinsen committed
                if not h1:
Yngve Levinsen's avatar
Yngve Levinsen committed
                    h1 = [thisdata[0] + " (std in paranthesis)"]
                    h2 = thisdata[2:10]
Yngve Levinsen's avatar
Yngve Levinsen committed
                d1.append(thisdata[1].split())
                d2.append(thisdata[10])
                d3.append(thisdata[11])

            # fix d1:
            for i in range(len(d1)):
                for j in range(len(d1[0])):
Yngve Levinsen's avatar
Yngve Levinsen committed
                    d1[i][j] = float(d1[i][j])
            d1 = np.array(d1)
            means = d1.mean(axis=0)
            stds = d1.std(axis=0)
            d1 = []
Yngve Levinsen's avatar
Yngve Levinsen committed
                if stds[i] / means[i] < 1e-10:
                    stds[i] = 0.0
Yngve Levinsen's avatar
Yngve Levinsen committed
                # some small std are removed..
Yngve Levinsen's avatar
Yngve Levinsen committed
                if stds[i] / means[i] > 1e-8:
                    d1.append("%f(%f)" % (means[i], stds[i]))
Yngve Levinsen's avatar
Yngve Levinsen committed
                else:  # error is 0
Yngve Levinsen's avatar
Yngve Levinsen committed
                    d1.append(str(means[i]))
            d1 = [" ".join(d1)]
Yngve Levinsen's avatar
Yngve Levinsen committed

            # create data:
Yngve Levinsen's avatar
Yngve Levinsen committed
            data = h1 + d1 + h2 + d2 + d3
                open(filename, "w").write("\n".join(data))

class envDiag:
    """
    Read ENV_diag1.dat file

    This contains e.g. the absolute phase at each diag

    For now we do not read in all info from the file,
    so feel free to request or add anything else you would like.


    def __init__(self, filename):
        self.filename = filename
        self.elements = {}
        # Needed to get an ordered dictionary:
        self._elementList = []
        self._readAsciiFile()
        self.units = {}
        self._setUnits()

    def _readAsciiFile(self):
        for line in open(self.filename, "r"):
            if lsp[0] == "DIAG":
                self.elements[int(lsp[2])] = {}
                self._elementList.append(int(lsp[2]))
                current = self.elements[int(lsp[2])]
                current["loc"] = float(lsp[4])
            elif lsp[0] == "Ibeam:":
                current["current"] = float(lsp[1])
            elif lsp[0] == "Positions":
                current["phase"] = float(lsp[5])
                current["energy"] = float(lsp[6])
            elif lsp[0] == "RMS":
                current["x_rms"] = float(lsp[4]) * 0.01
                current["y_rms"] = float(lsp[5]) * 0.01
                current["phase_rms"] = float(lsp[6])
                current["energy_rms"] = float(lsp[7])
            elif lsp[0] == "Emittances":
                current["emit_x"] = float(lsp[3])
                current["emit_y"] = float(lsp[4])
                current["emit_z"] = float(lsp[5])
            elif lsp[0] == "Emittances99":
                current["emit99_x"] = float(lsp[3])
                current["emit99_y"] = float(lsp[4])
                current["emit99_z"] = float(lsp[5])
            elif lsp[0] == "Twiss":
                if lsp[1] == "Alpha" and lsp[3] == "(XXp,":
                    current["alpha_x"] = float(lsp[6])
                    current["alpha_y"] = float(lsp[7])
                elif lsp[1] == "Alpha" and lsp[3] == "(ZDp/p)":
                    current["alpha_z"] = float(lsp[5])
                elif lsp[1] == "Beta":
                    current["beta_x"] = float(lsp[5])
                    current["beta_y"] = float(lsp[6])
                    current["beta_z"] = float(lsp[7])
        Set the units for each element in the element dictionary
        (empty string for all undefined)
        """
        for key in ["loc", "x_rms", "y_rms"]:
            self.units[key] = "m"
        for key in ["emit_x", "emit_y", "emit_z", "emit99_x", "emit99_y", "emit99_z"]:
            self.units[key] = "Pi.mm.mrad"
        for key in ["current"]:
            self.units[key] = "mA"
        for key in ["energy", "energy_rms"]:
            self.units[key] = "MeV"
        for key in ["phase", "phase_rms"]:
            self.units[key] = "deg"
        for key in ["beta_x", "beta_y", "beta_z"]:
            self.units[key] = "mm/mrad"

        for element in self.elements:
            for key in self.elements[element]:
                if key not in self.units:
                    self.units[key] = ""
        Make a pretty print of the content
        first = True
        rjust = 12
        for ekey in self._elementList:
            element = self.elements[ekey]

            # Print header if this is the first element..
            if first:
                keys = [key for key in element]
                print("#", end=" ")
                print("NUM".rjust(rjust), end=" ")
                    print(key.rjust(rjust), end=" ")
                print("#", end=" ")
                print("".rjust(rjust), end=" ")
                    print(self.units[key].rjust(rjust), end=" ")
            print("  " + str(ekey).rjust(rjust), end=" ")
            for key in keys:
                num = element[key]
                if isinstance(num, float):
                    strnum = "{:.5e}".format(num)
                print(strnum.rjust(rjust), end=" ")
            print()

    def getElement(self, elementId):
        Returns the element dictionary for the given ID
        return self.elements[elementId]

    def getElementAtLoc(self, location):
        Returns a list of elements at the location requested
        ret = []
        for key in self.elements:
            if abs(self.elements[key]["loc"] - location) < 1e-6:
                ret.append(self.elements[key])
        return ret

    def getParameterFromAll(self, parameter):
        Returns a list containing the given parameter from all DIAGS,
        ordered by the location of the DIAGs
        """
        if parameter == "NUM":
        ret = []
        for key in self._elementList:
            ret.append(self.elements[key][parameter])

        return ret
class partran(dict):
    Read partran1.out files..

    This class can also read tracewin.out (same format)
Yngve Levinsen's avatar
Yngve Levinsen committed
    def __init__(self, filename):
        self.filename = filename
        self._readAsciiFile()

    def _readAsciiFile(self):

        import numpy

        stream = open(self.filename, "r")
Yngve Levinsen's avatar
Yngve Levinsen committed
        for i in range(10):
Yngve Levinsen's avatar
Yngve Levinsen committed
            line = stream.readline()
            if line.strip()[0] == "#":
Yngve Levinsen's avatar
Yngve Levinsen committed
                break
        self.columns = ["NUM"] + line.split()[1:]
Yngve Levinsen's avatar
Yngve Levinsen committed
        self.data = numpy.loadtxt(stream)
Yngve Levinsen's avatar
Yngve Levinsen committed
        self._dict = {}
Yngve Levinsen's avatar
Yngve Levinsen committed
        for i in range(len(self.columns)):
Yngve Levinsen's avatar
Yngve Levinsen committed
            self[self.columns[i]] = self.data[:, i]

    Class to read in the field map structures

    WARNING: Work in progress!!

    def __init__(self, filename):
Yngve Levinsen's avatar
Yngve Levinsen committed
        self._filename = filename
        self._load_data(filename)

Yngve Levinsen's avatar
Yngve Levinsen committed
    def _load_data(self, filename):
        import os
        import numpy

        if not os.path.isfile(filename):
            raise ValueError("Cannot find file {}".format(filename))
        fin = open(filename, "r")
Yngve Levinsen's avatar
Yngve Levinsen committed
        line = fin.readline().split()
        self.header = []
Yngve Levinsen's avatar
Yngve Levinsen committed
        self.start = []
        self.end = []
        numindexes = []
Yngve Levinsen's avatar
Yngve Levinsen committed
        while len(line) > 1:
            [self.header.append(float(i)) for i in line]
            numindexes.append(int(line[0]) + 1)
            if len(line) == 2:
                self.start.append(0.0)
Yngve Levinsen's avatar
Yngve Levinsen committed
                self.end.append(float(line[1]))
Yngve Levinsen's avatar
Yngve Levinsen committed
                self.start.append(float(line[1]))
                self.end.append(float(line[2]))
            line = fin.readline().split()
        if len(self.start) == 1:
            self.z = numpy.mgrid[self.start[0] : self.end[0] : numindexes[0] * 1j]
Yngve Levinsen's avatar
Yngve Levinsen committed
            print(self.z)
        elif len(self.start) == 2:
            self.z, self.x = numpy.mgrid[
Yngve Levinsen's avatar
Yngve Levinsen committed
                self.start[0] : self.end[0] : numindexes[0] * 1j, self.start[1] : self.end[1] : numindexes[1] * 1j,
        elif len(self.start) == 3:
            self.z, self.x, self.y = numpy.mgrid[
                self.start[0] : self.end[0] : numindexes[0] * 1j,
                self.start[1] : self.end[1] : numindexes[1] * 1j,
                self.start[2] : self.end[2] : numindexes[2] * 1j,
            ]
Yngve Levinsen's avatar
Yngve Levinsen committed
        self.norm = float(line[0])
        self.header.append(self.norm)
Yngve Levinsen's avatar
Yngve Levinsen committed
        self.map = numpy.loadtxt(fin).reshape(numindexes)
    def get_flat_fieldmap(self):
        totmapshape = 1
        for i in self.map.shape:
            totmapshape *= i
        return self.map.reshape(totmapshape)

    def interpolate(self, npoints: tuple, method="cubic"):
        """
        Interpolate the map into a new mesh
        Each value should be an integer with the number of mesh points in each dimension
        intervals should be tuple-like with same number of elements
        as the map dimension, e.g. [0.8,0.8] for 2D
        Can also be a float if you want same interpolation factor in all planes

        method can be 'linear', 'nearest' or 'cubic'
        import numpy
        from scipy.interpolate import griddata

        values = self.map.flatten()

        if len(self.start) == 1:
            points = self.z[:]
            self.z = numpy.mgrid[self.start[0] : self.end[0] : npoints[0] * 1j]
            self.map = griddata(points, values, self.z)
        if len(self.start) == 2:
            points = numpy.array([self.z.flatten(), self.x.flatten()]).transpose()
            self.z, self.x = numpy.mgrid[
Yngve Levinsen's avatar
Yngve Levinsen committed
                self.start[0] : self.end[0] : npoints[0] * 1j, self.start[1] : self.end[1] : npoints[1] * 1j,
            ]
            self.map = griddata(points, values, (self.z, self.x))
        if len(self.start) == 3:
Yngve Levinsen's avatar
Yngve Levinsen committed
            points = numpy.array([self.z.flatten(), self.x.flatten(), self.y.flatten()]).transpose()
            self.z, self.x, self.y = numpy.mgrid[
                self.start[0] : self.end[0] : npoints[0] * 1j,
                self.start[1] : self.end[1] : npoints[1] * 1j,
                self.start[2] : self.end[2] : npoints[2] * 1j,
            ]
            self.map = griddata(points, values, (self.z, self.x, self.y))
            self.header[0] = npoints[0] - 1
            self.header[2] = npoints[1] - 1
            self.header[5] = npoints[2] - 1
    def savemap(self, filename):
        fout = open(filename, "w")
Yngve Levinsen's avatar
Yngve Levinsen committed
        for n, s in zip(self.map.shape, self.size):
            fout.write("{} {}\n".format(n - 1, s))
        fout.write("{}\n".format(self.norm))
Yngve Levinsen's avatar
Yngve Levinsen committed
        totmapshape = 1
        for i in self.map.shape:
Yngve Levinsen's avatar
Yngve Levinsen committed
            totmapshape *= i
        data = self.map.reshape(totmapshape)
            fout.write("{}\n".format(j))