dataMerge

#!/usr/bin/env python
# -*- coding: utf-8 -*-
"""
This program is part of pygimli
See https://www.pygimli.org for further information.
"""

import sys
import numpy as np

try:
    import pygimli as pg
    from pygimli.utils import trimDocString
except ImportError:
    sys.stderr.write('ERROR: cannot import the library pygimli.' +
                     ' Ensure that pygimli is in your PYTHONPATH')
    sys.exit(1)


def merge(data1, data2, ContainerType=pg.DataContainer, snap=0.001):
    """
        Merge two datacontainers into one by copying the sensor positions
        and datapoints from data2 into data1.\n
        Double sensor positions will be unified and snapped to a grid.
    """
    data = ContainerType(data1)
    data.add(data2, snap)
    return data


def loadProjectFile(
        projectfile, ContainerType=pg.DataContainer, verbose=False):
    """
        A project file defines how multiple data files can be imported and
        merged.
        A list of row entries for each data file to merge. The data will not be
        sorted, i.e., the sensor positions will be in the same order than the
        list entries are. However, sensors on the same positions (with
        tolerance snapping) will be merged.
        You can also give just one data file if you want just a spatial
        transformation of the data.

        The currently supported formats are:

        dataFilename
        dataFilename interpolationFilename
        dataFilename startx endx
        dataFilename startx starty endx endy
        dataFilename x1 y1 x2 y2 x3 y3 ...

        You can comment out a row by adding a leading '#'

        interpolationFilename is a 3-column-ascii-file (dx x y)

    """
    dataList = []

    fi = open(projectfile, "r")
    content = fi.readlines()
    fi.close()

    for c in content:
        row = c.split('\n')[0].split('#')[0].split()
        d = None

        if len(row) > 0 and row[0] != '#':
            if len(row) == 1:  # filename only
                d = ContainerType(row[0])
            elif len(row) == 2:  # dataFilename interpolationFilename
                # Thomas?? ist das von dir?? was macht das
                # kommt mir nur vage bekannt vor, benutzt aber extra-File
                d = ContainerType(row[0])

                xn = pg.x(d.sensorPositions())
                zn = pg.z(d.sensorPositions())

                tape, xt, yt = np.loadtxt(row[1], unpack=True)

                x3n = np.interp(xn, tape, xt)
                y3n = np.interp(xn, tape, yt)

                for i in range(d.sensorCount()):
                    d.setSensorPosition(i, pg.RVector3(x3n[i], y3n[i], zn[i]))

            elif len(row) == 3:  # filename xstart xend
                d = ContainerType(row[0])

                start = pg.RVector3(float(row[1]), 0.0)
                end = pg.RVector3(float(row[2]), 0.0)

                for i in range(d.sensorCount()):
                    pos = start + \
                        float(i) * (end - start) / (d.sensorCount() - 1.)
                    d.setSensorPosition(i, pos)

            elif len(row) == 5:  # filename xstart ystart xend yend
                d = ContainerType(row[0])

                start = pg.RVector3(float(row[1]), float(row[2]))
                end = pg.RVector3(float(row[3]), float(row[4]))

                for i in range(d.sensorCount()):
                    pos = start + \
                        float(i) * (end - start) / (d.sensorCount() - 1.)
                    d.setSensorPosition(i, pos)
            elif not len(row) % 2 is 0:  # dataFilename x1 y1 x2 y2 x3 y3 ...
                d = ContainerType(row[0])

                raise ('toimpl')
                start = pg.RVector3(
                    float(
                        row[1]), float(
                        row[2]), float(
                        row[3]))
                end = pg.RVector3(float(row[4]), float(row[5]), float(row[6]))

                for i in range(d.sensorCount()):
                    pos = start + \
                        float(i) * (end - start) / (d.sensorCount() - 1.)
                    d.setSensorPosition(i, pos)
            else:
                print(("cannot interpret the project format: len(row) = ",
                       len(row)))
                return dataList

            dataList.append(d)

            if verbose:
                print(("append: ", d))
                print(("from:", d.sensorPositions()[0],
                       "to:", d.sensorPositions()[-1]))

    return dataList


def main(argv):

    import argparse

    parser = argparse.ArgumentParser(
        formatter_class=argparse.RawDescriptionHelpFormatter,
        description=trimDocString(loadProjectFile.__doc__))

    parser.add_argument("-v", "--verbose", dest="verbose", action="store_true",
                        help="Be verbose.")
    parser.add_argument("-o", "--output", dest="outFileName", metavar="File",
                        help="Filename for the resulting data file.")
    parser.add_argument("-s", "--snap", dest="snap", type=float,
                        default=0.1,
                        help="Snap coordinates to gridsize")
    parser.add_argument("-B", "--bert", dest="bert", action="store_true",
                        help="Use BERT sensor indices (a b m n)")
    parser.add_argument('project')

    args = parser.parse_args()

    projectFileName = args.project

    if args.outFileName is None:
        args.outFileName = \
            projectFileName[0:projectFileName.find('.pro')] + '.dat'

    if args.verbose:
        print(args)
        print("verbose =", args.verbose)
        print("project =", projectFileName)
        print("output =", args.outFileName)
        print("snap =", args.snap)

    ContainerTyp = pg.DataContainer

    if args.bert:
        import pybert as b
        ContainerTyp = b.DataContainerERT

    dataList = loadProjectFile(projectFileName,
                               ContainerTyp,
                               verbose=args.verbose)

    outdata = dataList[0]

    if args.verbose:
        print("start merging ...")
        print(outdata)

    for d in dataList[1:]:
        outdata = merge(outdata, d, ContainerTyp, args.snap)

        if args.verbose:
            print(outdata)

    if args.verbose:
        print("Write file to: ", args.outFileName)

    outdata.save(args.outFileName)

if __name__ == "__main__":
    main(sys.argv[1:])