Commit f911ed1e authored by Germain Clavier's avatar Germain Clavier
Browse files

It's finally over! All files are python3 compatible (still need testing for debugging...)

parent 3f727e6e
Loading
Loading
Loading
Loading
+328 −303
Original line number Diff line number Diff line
@@ -5,12 +5,29 @@
# DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
# certain rights in this software.  This software is distributed under
# the GNU General Public License.

#
# matlab tool
#
#
#
# History
#   8/05, Matt Jones (BYU): original version
#
# ToDo list
#   allow choice of JPG or PNG or GIF when saving via "saveas" command
#   have vec/array import/export functions to pass variables between
#     Pizza.py and MatLab and have them named in MatLab

oneline = "Create plots via MatLab numerical analysis program"
# Variables
#   current = index of current figure (1-N)
#   figures = list of figure objects with each plot's attributes
#             so they aren't lost between replots
#   import command to yank MatLab variables back to Python
#   allow for alternate export command to name the variables from Python
#     in MatLab and vice versa

docstr = """
"""
Create plots via MatLab numerical analysis program
m = matlab()               start up MatLab
m.stop()                   shut down MatLab process

@@ -75,37 +92,24 @@ m.curve(N,'b','-','v') set color, line style, symbol of curve N
           '<' = left triangle, 'p' = pentagram, 'h' = hexagram
"""

# History
#   8/05, Matt Jones (BYU): original version

# ToDo list
#   allow choice of JPG or PNG or GIF when saving via "saveas" command
#   have vec/array import/export functions to pass variables between
#     Pizza.py and MatLab and have them named in MatLab

# Variables
#   current = index of current figure (1-N)
#   figures = list of figure objects with each plot's attributes
#             so they aren't lost between replots
#   import command to yank MatLab variables back to Python
#   allow for alternate export command to name the variables from Python
#     in MatLab and vice versa

# Imports and external programs

import types, os
import os

try: from DEFAULTS import PIZZA_MATLAB
except: PIZZA_MATLAB = "matlab -nosplash -nodesktop -nojvm"
try:
    from DEFAULTS import PIZZA_MATLAB
except (ModuleNotFoundError, ImportError):
    PIZZA_MATLAB = "matlab -nosplash -nodesktop -nojvm"

# Class definition


class matlab:

    # --------------------------------------------------------------------

    def __init__(self):
    self.MATLAB = os.popen(PIZZA_MATLAB,'w')
        self.MATLAB = os.popen(PIZZA_MATLAB, "w")
        self.file = "tmp.matlab"
        self.figures = []
        self.select(1)
@@ -119,15 +123,16 @@ class matlab:
    # --------------------------------------------------------------------

    def __call__(self, command):
    self.MATLAB.write(command + '\n')
        self.MATLAB.write(command + "\n")
        self.MATLAB.flush()

    # --------------------------------------------------------------------

    def enter(self):
        while 1:
      command = raw_input("matlab> ")
      if command == "quit" or command == "exit": return
            command = input("matlab> ")
            if command == "quit" or command == "exit":
                return
            self.__call__(command)

    # --------------------------------------------------------------------
@@ -136,11 +141,12 @@ class matlab:
    def plot(self, *vectors):
        if len(vectors) == 1:
            file = self.file + ".%d.1" % self.current
      linear = range(len(vectors[0]))
            linear = list(range(len(vectors[0])))
            self.export(file, linear, vectors[0])
            self.figures[self.current - 1].ncurves = 1
        else:
      if len(vectors) % 2: raise StandardError,"vectors must come in pairs"
            if len(vectors) % 2:
                sys.exit("vectors must come in pairs")
            for i in range(0, len(vectors), 2):
                file = self.file + ".%d.%d" % (self.current, i / 2 + 1)
                self.export(file, vectors[i], vectors[i + 1])
@@ -156,14 +162,20 @@ class matlab:
        for i in range(start, stop, skip):
            partial_vecs = []
            for vec in vectors:
        if i: partial_vecs.append(vec[:i])
        else: partial_vecs.append([0])
                if i:
                    partial_vecs.append(vec[:i])
                else:
                    partial_vecs.append([0])
            self.plot(*partial_vecs)

      if n < 10:     newfile = file + "000" + str(n)
      elif n < 100:  newfile = file + "00" + str(n)
      elif n < 1000: newfile = file + "0" + str(n)
      else:          newfile = file + str(n)
            if n < 10:
                newfile = file + "000" + str(n)
            elif n < 100:
                newfile = file + "00" + str(n)
            elif n < 1000:
                newfile = file + "0" + str(n)
            else:
                newfile = file + str(n)

            self.save(newfile)
            n += 1
@@ -174,13 +186,14 @@ class matlab:
    def export(self, filename, *vectors):
        n = len(vectors[0])
        for vector in vectors:
      if len(vector) != n: raise StandardError,"vectors must be same length"
    f = open(filename,'w')
            if len(vector) != n:
                sys.exit("vectors must be same length")
        f = open(filename, "w")
        nvec = len(vectors)
    for i in xrange(n):
      for j in xrange(nvec):
        print >>f,vectors[j][i],
      print >>f
        for i in range(n):
            for j in range(nvec):
                f.write("{}".format(vectors[j][i]))
            f.write('\n')
        f.close()

    # --------------------------------------------------------------------
@@ -207,11 +220,13 @@ class matlab:
    #   use tmp.done as semaphore to indicate plot is finished

    def save(self, file):
    if os.path.exists("tmp.done"): os.remove("tmp.done")
        if os.path.exists("tmp.done"):
            os.remove("tmp.done")
        cmd = "saveas(gcf,'%s.eps','psc2')" % file
        self.__call__(cmd)
        self.__call__("!touch tmp.done")
    while not os.path.exists("tmp.done"): continue
        while not os.path.exists("tmp.done"):
            continue

    # --------------------------------------------------------------------
    # restore default attributes by creating a new fig object
@@ -305,10 +320,13 @@ class matlab:

    def curve(self, num, *settings):
        fig = self.figures[self.current - 1]
    while len(fig.curves) < num: fig.curves.append('')
        while len(fig.curves) < num:
            fig.curves.append("")
        value = settings[0]
    if len(settings) >= 2: value += settings[1]
    if len(settings) >= 3: value += settings[2]
        if len(settings) >= 2:
            value += settings[1]
        if len(settings) >= 3:
            value += settings[2]
        fig.curves[num - 1] = value
        self.draw()

@@ -318,15 +336,20 @@ class matlab:

    def draw(self):
        fig = self.figures[self.current - 1]
    if not fig.ncurves: return
        if not fig.ncurves:
            return
        self.__call__("figure(%d)" % self.current)

    if not fig.xlog and not fig.ylog: cmd = "plot"
    elif fig.xlog and not fig.ylog: cmd = "semilogx"
    elif not fig.xlog and fig.ylog: cmd = "semilogy"
    elif fig.xlog and fig.ylog: cmd = "loglog"
        if not fig.xlog and not fig.ylog:
            cmd = "plot"
        elif fig.xlog and not fig.ylog:
            cmd = "semilogx"
        elif not fig.xlog and fig.ylog:
            cmd = "semilogy"
        elif fig.xlog and fig.ylog:
            cmd = "loglog"

    cmd += '('
        cmd += "("
        for i in range(fig.ncurves):
            file = self.file + ".%d.%d" % (self.current, i + 1)
            readcmd = "pizza%d = importdata('%s');" % (i + 1, file)
@@ -347,7 +370,8 @@ class matlab:
        self.__call__("xlabel('%s','FontSize',16)" % fig.xtitle)
        self.__call__("ylabel('%s','FontSize',16)" % fig.ytitle)

    if fig.xlimit == 0 or fig.ylimit == 0: self.__call__("axis auto")
        if fig.xlimit == 0 or fig.ylimit == 0:
            self.__call__("axis auto")
        if fig.xlimit:
            self.__call__("xlim([%g,%g])" % (fig.xlimit[0], fig.xlimit[1]))
        if fig.ylimit:
@@ -359,11 +383,12 @@ class matlab:
            text = fig.labels[i][2]  # kludge to set label font size
            self.__call__("text(%g,%g,'%s','FontSize',16)" % (x, y, text))


# --------------------------------------------------------------------
# class to store settings for a single plot

class figure:

class figure:
    def __init__(self):
        self.ncurves = 0
        self.curves = []
+797 −541

File changed.

Preview size limit exceeded, changes collapsed.

+1643 −1404

File changed.

Preview size limit exceeded, changes collapsed.

+752 −627

File changed.

Preview size limit exceeded, changes collapsed.

+344 −323
Original line number Diff line number Diff line
@@ -5,12 +5,31 @@
# DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
# certain rights in this software.  This software is distributed under
# the GNU General Public License.

#
# tdump tool
#
#
#
# History
#   4/11, Steve Plimpton (SNL): original version
#
# Variables
#   flist = list of dump file names
#   increment = 1 if reading snapshots one-at-a-time
#   nextfile = which file to read from via next()
#   eof = ptr into current file for where to read via next()
#   nsnaps = # of snapshots
#   snaps = list of snapshots
#   names = dictionary of column names:
#     key = "id", value = column # (0 to M-1)
#   Snap = one snapshot
#     time = time stamp
#     natoms = # of atoms
#     xlo,xhi,ylo,yhi,zlo,zhi = box bounds (float)
#     atoms[i][j] = 2d array of floats, i = 0 to natoms-1, j = 0 to ncols-1

oneline = "Read dump files with triangle info"

docstr = """
"""
Read dump files with triangle info
t = tdump("dump.one")              read in one or more dump files
t = tdump("dump.1 dump.2.gz")      can be gzipped
t = tdump("dump.*")                wildcard expands to multiple files
@@ -36,7 +55,7 @@ time,box,atoms,bonds,tris,lines = t.viz(index) return list of viz objects
  viz() returns line info for specified timestep index
    can also call as viz(time,1) and will find index of preceding snapshot
    time = timestep value
    box = \[xlo,ylo,zlo,xhi,yhi,zhi\]
    box = [xlo,ylo,zlo,xhi,yhi,zhi]
    atoms = NULL
    bonds = NULL
    tris = id,type,x1,y1,z1,x2,y2,z2,x3,y3,z3 for each tri as 2d array
@@ -49,60 +68,44 @@ t.owrap(...) wrap tris to same image as their atoms
  useful for wrapping all molecule's atoms/tris the same so it is contiguous
"""

# History
#   4/11, Steve Plimpton (SNL): original version

# Variables
#   flist = list of dump file names
#   increment = 1 if reading snapshots one-at-a-time
#   nextfile = which file to read from via next()
#   eof = ptr into current file for where to read via next()
#   nsnaps = # of snapshots
#   snaps = list of snapshots
#   names = dictionary of column names:
#     key = "id", value = column # (0 to M-1)
#   Snap = one snapshot
#     time = time stamp
#     natoms = # of atoms
#     xlo,xhi,ylo,yhi,zlo,zhi = box bounds (float)
#     atoms[i][j] = 2d array of floats, i = 0 to natoms-1, j = 0 to ncols-1

# Imports and external programs

import sys, commands, re, glob, types
import sys
import subprocess
import re
import glob
import types
from math import sqrt
from os import popen

try:
import numpy as np
    oldnumeric = False
except:
    import Numeric as np
    oldnumeric = True

try: from DEFAULTS import PIZZA_GUNZIP
except: PIZZA_GUNZIP = "gunzip"
try:
    from DEFAULTS import PIZZA_GUNZIP
except (ModuleNotFoundError, ImportError):
    PIZZA_GUNZIP = "gunzip"

# Class definition


class tdump:

    # --------------------------------------------------------------------

  def __init__(self,*list):
    def __init__(self, *args):
        self.snaps = []
        self.nsnaps = 0
        self.names = {}

        # flist = list of all dump file names

    words = list[0].split()
        words = args[0].split()
        self.flist = []
    for word in words: self.flist += glob.glob(word)
    if len(self.flist) == 0 and len(list) == 1:
      raise StandardError,"no ldump file specified"
        for word in words:
            self.flist += glob.glob(word)
        if len(self.flist) == 0 and len(args) == 1:
            sys.exit("no ldump file specified")

    if len(list) == 1:
        if len(args) == 1:
            self.increment = 0
            self.read_all()
        else:
@@ -119,44 +122,47 @@ class tdump:

        for file in self.flist:
            if file[-3:] == ".gz":
        f = popen("%s -c %s" % (PIZZA_GUNZIP,file),'r')
      else: f = open(file)
                f = popen("%s -c %s" % (PIZZA_GUNZIP, file), "r")
            else:
                f = open(file)

            snap = self.read_snapshot(f)
            while snap:
                self.snaps.append(snap)
        print snap.time,
                print(snap.time)
                sys.stdout.flush()
                snap = self.read_snapshot(f)

            f.close()
    print
        print("\n")

        # sort entries by timestep, cull duplicates

        self.snaps.sort(self.compare_time)
        self.cull()
        self.nsnaps = len(self.snaps)
    print "read %d snapshots" % self.nsnaps
        print("read %d snapshots" % self.nsnaps)

    # --------------------------------------------------------------------
    # read next snapshot from list of files

    def next(self):

    if not self.increment: raise StandardError,"cannot read incrementally"
        if not self.increment:
            sys.exit("cannot read incrementally")

        # read next snapshot in current file using eof as pointer
        # if fail, try next file
        # if new snapshot time stamp already exists, read next snapshot

        while 1:
      f = open(self.flist[self.nextfile],'rb')
            f = open(self.flist[self.nextfile], "rb")
            f.seek(self.eof)
            snap = self.read_snapshot(f)
            if not snap:
                self.nextfile += 1
	if self.nextfile == len(self.flist): return -1
                if self.nextfile == len(self.flist):
                    return -1
                f.close()
                self.eof = 0
                continue
@@ -165,7 +171,8 @@ class tdump:
            try:
                self.findtime(snap.time)
                continue
      except: break
            except Exception:
                break

        self.snaps.append(snap)
        snap = self.snaps[self.nsnaps]
@@ -198,21 +205,21 @@ class tdump:
            if snap.natoms:
                words = f.readline().split()
                ncol = len(words)
        for i in xrange(1,snap.natoms):
                for i in range(1, snap.natoms):
                    words += f.readline().split()
        floats = map(float,words)
        if oldnumeric: atoms = np.zeros((snap.natoms,ncol),np.Float)
        else: atoms = np.zeros((snap.natoms,ncol),np.float)
                values = list(map(float, words))
                atoms = np.zeros((snap.natoms, ncol), dtype=np.float)
                start = 0
                stop = ncol
        for i in xrange(snap.natoms):
          atoms[i] = floats[start:stop]
                for i in range(snap.natoms):
                    atoms[i] = values[start:stop]
                    start = stop
                    stop += ncol
      else: atoms = None
            else:
                atoms = None
            snap.atoms = atoms
            return snap
    except:
        except Exception:
            return 0

    # --------------------------------------------------------------------
@@ -220,7 +227,7 @@ class tdump:

    def map(self, *pairs):
        if len(pairs) % 2 != 0:
      raise StandardError, "tdump map() requires pairs of mappings"
            sys.exit("tdump map() requires pairs of mappings")
        for i in range(0, len(pairs), 2):
            j = i + 1
            self.names[pairs[j]] = pairs[i] - 1
@@ -250,9 +257,10 @@ class tdump:
    # --------------------------------------------------------------------

    def findtime(self, n):
    for i in xrange(self.nsnaps):
      if self.snaps[i].time == n: return i
    raise StandardError, "no step %d exists" % n
        for i in range(self.nsnaps):
            if self.snaps[i].time == n:
                return i
        sys.exit("no step %d exists" % n)

    # --------------------------------------------------------------------
    # delete successive snapshots with duplicate time stamp
@@ -270,21 +278,23 @@ class tdump:
    # if called with flag, then index is timestep, so convert to snapshot index

    def viz(self, index, flag=0):
    if not flag: isnap = index
        if not flag:
            isnap = index
        else:
            times = self.time()
            n = len(times)
            i = 0
            while i < n:
        if times[i] > index: break
                if times[i] > index:
                    break
                i += 1
            isnap = i - 1
        snap = self.snaps[isnap]

        time = snap.time
        box = [snap.xlo, snap.ylo, snap.zlo, snap.xhi, snap.yhi, snap.zhi]
    id = self.names["id"]
    type = self.names["type"]
        atid = self.names["id"]
        attype = self.names["type"]
        corner1x = self.names["corner1x"]
        corner1y = self.names["corner1y"]
        corner1z = self.names["corner1z"]
@@ -299,15 +309,22 @@ class tdump:
        # don't add tri if all 4 values are 0 since not a line

        tris = []
    for i in xrange(snap.natoms):
        for i in range(snap.natoms):
            atom = snap.atoms[i]
            c1 = [atom[corner1x], atom[corner1y], atom[corner1z]]
            c2 = [atom[corner2x], atom[corner2y], atom[corner2z]]
            c3 = [atom[corner3x], atom[corner3y], atom[corner3z]]
            n = normal(c1, c2, c3)
      if c1[0] == 0.0 and c1[1] == 0.0 and c1[2] == 0.0 and \
              c2[0] == 0.0 and c2[1] == 0.0 and c2[2] == 0.0: continue
      tris.append([atom[id],atom[type]] + c1 + c2 + c3 + n)
            if (
                c1[0] == 0.0
                and c1[1] == 0.0
                and c1[2] == 0.0
                and c2[0] == 0.0
                and c2[1] == 0.0
                and c2[2] == 0.0
            ):
                continue
            tris.append([atom[atid], atom[attype]] + c1 + c2 + c3 + n)

        return time, box, None, None, tris, None

@@ -316,7 +333,7 @@ class tdump:
    # invoked by dump() when it does an owrap() on its atoms

    def owrap(self, time, xprd, yprd, zprd, idsdump, atomsdump, iother, ix, iy, iz):
    id = self.names["id"]
        atid = self.names["id"]
        corner1x = self.names["corner1x"]
        corner1y = self.names["corner1y"]
        corner1z = self.names["corner1z"]
@@ -335,8 +352,8 @@ class tdump:
        # jdump = atom J in dump's atoms that atom I was owrapped on
        # delx,dely = offset applied to atom I and thus to line I

    for i in xrange(snap.natoms):
      tag = atoms[i][id]
        for i in range(snap.natoms):
            tag = atoms[i][atid]
            idump = idsdump[tag]
            jdump = idsdump[atomsdump[idump][iother]]
            delx = (atomsdump[idump][ix] - atomsdump[jdump][ix]) * xprd
@@ -352,15 +369,19 @@ class tdump:
            atoms[i][corner3y] += dely
            atoms[i][corner3z] += delz


# --------------------------------------------------------------------
# one snapshot


class Snap:
    pass


# --------------------------------------------------------------------
# compute normal for a triangle with 3 vertices


def normal(x, y, z):
    v1 = 3 * [0]
    v1[0] = y[0] - x[0]
Loading