program to resize simulations

py/mains/compare_snapshots.py

+#!/usr/bin/env python
+
+#
+# This program test unsio library by readind and saving same file in different output
+# file format (gadget2, nemo) and comparing all arrays with original one
+#
+from __future__ import print_function    
+import sys
+import argparse
+import numpy as np                # arrays are treated as numpy arrays
+
+import os.path
+#dirname, filename = os.path.split(os.path.abspath(__file__))
+#sys.path.append(dirname+'../modules/')  # trick to find modules directory
+from py_unsio import *
+import copy
+import tempfile
+#from IPython import embed
+
+class snap:
+    time   = None
+    nbody  = None
+    mass   = None
+    pos    = None
+    vel    = None
+    acc    = None
+    pot    = None
+    id     = None
+    age    = None
+    hsml   = None
+    rho    = None
+    metal  = None
+    interface = ""
+ 
+# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+# commandLine, parse the command line 
+def commandLine():
+    # help
+    parser = argparse.ArgumentParser(description="Compare to given snapshot arrays vs arrays",
+                                    formatter_class=argparse.ArgumentDefaultsHelpFormatter)
+    # options
+    parser.add_argument('snapshot_1', help='UNS snapshot')
+    parser.add_argument('snapshot_2', help='UNS snapshot')
+    parser.add_argument('component', help='component name')
+    float_parser=parser.add_mutually_exclusive_group(required=True)
+    float_parser.add_argument('--float',    help='floating format operation',action='store_true')
+    float_parser.add_argument('--double', help='double format operation',  action='store_true')
+    parser.add_argument('--out', help='save output file name ?',default="")
+        # parse
+    args = parser.parse_args()
+
+    # start main funciton
+    process(args)
+
+# -----------------------------------------------------
+def readSnap(simname, comp, single):
+    components=comp 
+    verbose=False
+
+    
+    #timef=float(times)
+    
+    # Create a UNSIO object
+    if (single) :
+        #print("Float format object")
+        uns = CunsIn(simname,components,"all",verbose)
+    else:
+        #print("Double format object")
+        uns = CunsInD(simname,components,"all",verbose)
+
+    #print ("simname=",simname,file=sys.stderr)
+
+    mysnap=snap() # instantiate a snap object
+    # load frame
+    ok=uns.nextFrame("")
+    #print ok
+
+    if (ok) :
+        #embed()
+        mysnap.interface = uns.getInterfaceType()
+        ok,mysnap.time = uns.getValueF("time")
+        ok,mysnap.pos  = uns.getArrayF(comp,"pos")
+        ok,mysnap.vel  = uns.getArrayF(comp,"vel")
+        ok,mysnap.mass = uns.getArrayF(comp,"mass")
+        ok,mysnap.hsml = uns.getArrayF(comp,"hsml")
+        ok,mysnap.rho  = uns.getArrayF(comp,"rho")
+        ok,mysnap.age  = uns.getArrayF(comp,"age")
+        ok,mysnap.acc  = uns.getArrayF(comp,"acc")
+        ok,mysnap.pot  = uns.getArrayF(comp,"pot")
+        ok,mysnap.metal= uns.getArrayF(comp,"metal")
+        if ok:
+            mysnap.metal=fixMetal(mysnap.metal)
+        ok,mysnap.id   = uns.getArrayI(comp,"id")
+        uns.close()
+        return True,copy.deepcopy(mysnap)
+        
+    else :
+        print ("Didn't load anything....",file=sys.stderr)
+
+    return False
+  
+# -----------------------------------------------------
+def compareArray(CA,CB,attr):
+    #embed()
+    A=getattr(CA,attr)
+    B=getattr(CB,attr)
+    ok=False
+    disp=True
+    if notCompare(CA,CB,attr):
+        disp=False
+    else:
+        if attr=="time":
+            ok = (A==B)
+        else:
+            ok=(A==B).all()
+            if ok :
+                if (A.size):
+                    disp=True
+                else:
+                    disp=False
+
+    if (disp):
+        print("[",attr,"]",ok)
+        if not ok:
+            print("\tA:",A[0:2],"\n\tB:",B[0:2])
+
+# -----------------------------------------------------
+# do not compare in the following cases
+def notCompare(CA,CB,attr):
+    status=False
+    if (CA.interface=="Nemo" or CB.interface=="Nemo"):
+        if attr=="metal":
+            status=True
+        if attr=="age":
+            status=True
+        if status:
+            print("<",attr,"> attribute not supported with NEMO format")
+    A=getattr(CA,attr)
+    if (attr != "time" and A.size==0):
+        status=True
+        #print("In <",attr,"> attribute not tested")
+    B=getattr(CB,attr)
+    if (attr != "time" and A.size!=0 and B.size==0):
+        print("In <",attr,"> attribute missing in B")
+        status=True
+    return status
+        
+# -----------------------------------------------------
+def compare(CA,CB,comp):
+    print("-----------------------------------------------------")
+    print("Comparing : <%s> [%s] vs [%s]\n"%(comp,CA.interface,CB.interface))
+    for attr in ("pos","vel","mass","age","hsml","rho","metal","acc","pot", "id","time"):
+        compareArray(CA,CB,attr)
+
+# -----------------------------------------------------
+def fixMetal(metal):
+    if (metal==-1.0).all() :
+        print("fixing metal....")
+        return np.empty(0)
+    else:
+        return metal
+
+# -----------------------------------------------------
+# process
+def process(args):
+
+    if args.component == "all":
+        comp_sel=["halo", "gas", "stars", "disk", "bulge", "bndry"]
+    else:
+        comp_sel=[args.component]
+    for comp in (comp_sel):
+        ok,insnap1=readSnap(args.snapshot_1,comp,args.float)
+        ok,insnap2=readSnap(args.snapshot_2,comp,args.float)
+        compare(insnap1,insnap2,comp)
+
+
+# -----------------------------------------------------
+# main program
+commandLine()   # parse command line
+#

py/mains/remove_sim_halo.py

+#!/usr/bin/env python
+from __future__ import print_function
+
+import sys
+sys.path=['/home/jcl/works/GIT/uns_projects/py/modules/','/home/jcl/works/GIT/uns_projects/py/modules/simulations']+sys.path
+
+#from py_unstools import *         # import py_unstools package
+from py_unsio import *
+from uns_simu import *
+
+
+import argparse
+import os,subprocess
+
+
+from simulations.creducesim import *
+
+# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+# commandLine, parse the command line 
+def commandLine():
+    dbname=None
+    modulo=10
+    test=False
+
+     # help
+    parser = argparse.ArgumentParser(description="Remove halo from simulation",
+                                    formatter_class=argparse.ArgumentDefaultsHelpFormatter)
+    # options
+    parser.add_argument('simname', help='UNS Simulation name')
+    parser.add_argument('--dir', help='directory to store new files',default=None)
+    parser.add_argument('--keep', help='keep halo every frequency',default=modulo,type=int)
+    parser.add_argument('--overwrite',help='overwrite new frame if present', default=False)
+    parser.add_argument('--test',help='test without doing anything', dest="test",action="store_true", default=test)
+    parser.add_argument('--dbname',help='UNS database file name', default=dbname)
+    parser.add_argument('--verbose',help='verbose mode',dest="verbose", action="store_true", default=False)
+     # parse
+    args = parser.parse_args()
+
+    # start main funciton
+    process(args)
+
+# -----------------------------------------------------
+# process, is the core function
+def process(args):
+    try:
+        simu=CReducesim(simname=args.simname,keep=args.keep,overwrite=args.overwrite,dir=args.dir,test=args.test,verbose=args.verbose)
+    except Exception as x :
+        print (x.message)
+    else:
+        simu.resizeSim()
+
+
+# -----------------------------------------------------
+# process, is the core function
+def process2(args):
+    try:
+        simu=UnsSimu(args.simname,verbose=args.verbose)
+    except Exception as x :
+        print (x.message)
+    else:
+        info=simu.getInfo() # get simulation info
+        simname=info["name"]
+        simtype=(info["dir"]).split("/")[2]
+        list=simu.getSnapshotList()
+
+        # get real path even if there is a link
+        # remove /net/direct in case mounted
+        realpath=(os.path.realpath(list[0])).replace("/net/direct","")
+        print("Real Path [%s]"%(realpath),file=sys.stderr)
+
+        if args.dir is None:
+            # get file system root
+            inrootdir=realpath.split("/")[1]
+            print("inrootdir = [%s] <%s>"%(inrootdir,simtype),file=sys.stderr)
+            args.dir="/"+inrootdir+"/"+simtype+"2/"+simname+"/SNAPS"
+            print("target dir <%s>"%(args.dir),file=sys.stderr)
+
+        # check dir
+        if not os.path.isdir(args.dir):
+            print("Create directory [%s]"%(args.dir),file=sys.stderr)
+            try:
+                if not args.test:
+                    os.makedirs(args.dir)
+            except :
+                print("Unable to create directory [%s]"%(args.dir))
+                sys.exit()
+        cpt=1
+        for snap in list:
+            print(snap,file=sys.stderr)
+            
+            if cpt%args.keep == 1:
+                select="all"
+            else:
+                select="gas,stars,disk,bndry,bulge"
+                pass
+            basename=os.path.basename(snap) # input
+            newsnap=args.dir+"/"+basename   # output
+            # check output file
+            if not os.path.isfile(newsnap) or args.overwrite:
+                cmd="uns2uns %s %s select=%s type=gadget3"%(snap,newsnap,select)
+                print("<%s>"%(cmd),file=sys.stderr)
+                if not args.test:
+                    subprocess.call(["uns2uns",snap,newsnap,"select="+select,"type=gadget3"])
+
+            cpt+=1
+
+# -----------------------------------------------------
+# main program
+if __name__ == '__main__':
+  commandLine()
+
+
+