223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

old shell scripts 

- removed unmaintained scripts from 'legacy'
- pushed deprecated scripts to 'legacy'
This commit is contained in:
Martin Diehl 2021-03-26 08:48:05 +01:00
parent 35c854b549
commit 75fb080638
21 changed files with 0 additions and 1273 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

55
.gitlab-ci.yml
View File

@ -105,39 +105,12 @@ pytest_python:
- release
###################################################################################################
Pre_SeedGeneration:
stage: deprecated
script: PreProcessing_SeedGeneration/test.py
except:
- master
- release
Pre_GeomGeneration:
stage: deprecated
script: PreProcessing_GeomGeneration/test.py
except:
- master
- release
Pre_GeomModification:
stage: deprecated
script: PreProcessing_GeomModification/test.py
except:
- master
- release
Pre_General:
stage: deprecated
script: PreProcessing/test.py
except:
- master
- release
Post_General:
stage: deprecated
script: PostProcessing/test.py
except:
- master
- release
Post_GeometryReconstruction:
@ -147,27 +120,6 @@ Post_GeometryReconstruction:
- master
- release
Post_addCurl:
stage: deprecated
script: addCurl/test.py
except:
- master
- release
Post_addDivergence:
stage: deprecated
script: addDivergence/test.py
except:
- master
- release
Post_addGradient:
stage: deprecated
script: addGradient/test.py
except:
- master
- release
Post_OrientationAverageMisorientation:
stage: deprecated
script:
@ -262,13 +214,6 @@ pytest_fortran:
- master
- release
Plasticity_DetectChanges:
stage: grid
script: Plasticity_DetectChanges/test.py
except:
- master
- release
Phenopowerlaw_singleSlip:
stage: grid
script: Phenopowerlaw_singleSlip/test.py

0
processing/post/DADF5toDREAM3D.py → processing/legacy/DADF5toDREAM3D.py
View File

51
processing/legacy/addAPS34IDEstrainCoords.py
View File

@ -1,51 +0,0 @@
#!/usr/bin/env python3
import os
import sys
from io import StringIO
from optparse import OptionParser
import numpy as np
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [ASCIItable(s)]', description = """
Transform X,Y,Z,F APS BeamLine 34 coordinates to x,y,z APS strain coordinates.
""", version = scriptID)
parser.add_option('-f','--frame',dest='frame', metavar='string',
help='label of APS X,Y,Z coords')
parser.add_option('--depth', dest='depth', metavar='string',
help='depth')
(options,filenames) = parser.parse_args()
if filenames == []: filenames = [None]
if options.frame is None:
parser.error('frame not specified')
if options.depth is None:
parser.error('depth not specified')
rot_to_TSL = damask.Rotation.from_axis_angle([-1,0,0,.75*np.pi])
for name in filenames:
damask.util.report(scriptName,name)
table = damask.Table.load(StringIO(''.join(sys.stdin.read())) if name is None else name)
coord = - table.get(options.frame)
coord[:,2] += table.get(options.depth)[:,0]
table.add('coord',rot_to_TSL.broadcast_to(coord.shape[0]) @ coord,scriptID+' '+' '.join(sys.argv[1:]))\
.save((sys.stdout if name is None else name),legacy=True)

0
processing/post/addCompatibilityMismatch.py → processing/legacy/addCompatibilityMismatch.py
View File

48
processing/legacy/addCumulative.py
View File

@ -1,48 +0,0 @@
#!/usr/bin/env python3
import os
import sys
from io import StringIO
from optparse import OptionParser
import numpy as np
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [ASCIItable(s)]', description = """
Add cumulative (sum of first to current row) values for given label(s).
""", version = scriptID)
parser.add_option('-l','--label',
dest='labels',
action = 'extend', metavar = '<string LIST>',
help = 'columns to cumulate')
parser.add_option('-p','--product',
dest='product', action = 'store_true',
help = 'product of values instead of sum')
(options,filenames) = parser.parse_args()
if filenames == []: filenames = [None]
if options.labels is None:
parser.error('no data column(s) specified.')
for name in filenames:
damask.util.report(scriptName,name)
table = damask.Table.load(StringIO(''.join(sys.stdin.read())) if name is None else name)
for label in options.labels:
table = table.add('cum_{}({})'.format('prod' if options.product else 'sum',label),
np.cumprod(table.get(label),0) if options.product else np.cumsum(table.get(label),0),
scriptID+' '+' '.join(sys.argv[1:]))
table.save((sys.stdout if name is None else name),legacy=True)

0
processing/post/addCurl.py → processing/legacy/addCurl.py
View File

0
processing/post/addDisplacement.py → processing/legacy/addDisplacement.py
View File

0
processing/post/addDivergence.py → processing/legacy/addDivergence.py
View File

100
processing/legacy/addEhkl.py
View File

@ -1,100 +0,0 @@
#!/usr/bin/env python3
import os
import sys
from optparse import OptionParser
import numpy as np
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
def normalize(vec):
return vec/np.sqrt(np.inner(vec,vec))
def E_hkl(stiffness,vec): # stiffness = (c11,c12,c44)
v = normalize(vec)
S11 = (stiffness[0]+stiffness[1])/(stiffness[0]*stiffness[0]+stiffness[0]*stiffness[1]-2.0*stiffness[1]*stiffness[1])
S12 = ( -stiffness[1])/(stiffness[0]*stiffness[0]+stiffness[0]*stiffness[1]-2.0*stiffness[1]*stiffness[1])
S44 = 1.0/stiffness[2]
invE = S11-(S11-S12-0.5*S44)* (1.0 - \
(v[0]**4+v[1]**4+v[2]**4) \
/#------------------------------------
np.inner(v,v)**2 \
)
return 1.0/invE
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [ASCIItable(s)]', description = """
Add column(s) containing directional stiffness based on given cubic stiffness values C11, C12, and C44 in consecutive columns.
""", version = scriptID)
parser.add_option('-c','--stiffness',
dest = 'stiffness',
action = 'extend', metavar = '<string LIST>',
help = 'heading of column containing C11 (followed by C12, C44) field values')
parser.add_option('-d','--direction','--hkl',
dest = 'hkl',
type = 'int', nargs = 3, metavar = 'int int int',
help = 'direction of elastic modulus [%default]')
parser.set_defaults(hkl = (1,1,1),
)
(options,filenames) = parser.parse_args()
if options.stiffness is None:
parser.error('no data column specified...')
# --- loop over input files -------------------------------------------------------------------------
if filenames == []: filenames = [None]
for name in filenames:
try:
table = damask.ASCIItable(name = name)
except IOError:
continue
damask.util.report(scriptName,name)
# ------------------------------------------ read header ------------------------------------------
table.head_read()
# ------------------------------------------ sanity checks ----------------------------------------
remarks = []
columns = []
for i,column in enumerate(table.label_index(options.stiffness)):
if column < 0: remarks.append('column {} not found.'.format(options.stiffness[i]))
else:
columns.append(column)
table.labels_append(['E{}{}{}({arg2})'.format(*options.hkl,arg2=options.stiffness[i])]) # extend ASCII header with new labels
if remarks != []: damask.util.croak(remarks)
# ------------------------------------------ assemble header --------------------------------------
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
table.head_write()
# ------------------------------------------ process data ------------------------------------------
outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table
for column in columns:
table.data_append(E_hkl(list(map(float,table.data[column:column+3])),options.hkl))
outputAlive = table.data_write() # output processed line
# ------------------------------------------ output finalization -----------------------------------
table.close() # close ASCII tables

0
processing/post/addEuclideanDistance.py → processing/legacy/addEuclideanDistance.py
View File

0
processing/post/addGradient.py → processing/legacy/addGradient.py
View File

0
processing/post/addSchmidfactors.py → processing/legacy/addSchmidfactors.py
View File

0
processing/pre/geom_grainGrowth.py → processing/legacy/geom_grainGrowth.py
View File

0
processing/pre/gmsh_identifySurfaces.py → processing/legacy/gmsh_identifySurfaces.py
View File

194
processing/legacy/imageDataDeformed.py
View File

@ -1,194 +0,0 @@
#!/usr/bin/env python3
import os
import sys
from optparse import OptionParser
import numpy as np
from PIL import Image, ImageDraw
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [file[s]]', description = """
Generate PNG image from scalar data on grid deformed by (periodic) deformation gradient.
""", version = scriptID)
parser.add_option('-l','--label',
dest = 'label',
type = 'string', metavar = 'string',
help = 'column containing data [all]')
parser.add_option('-r','--range',
dest = 'range',
type = 'float', nargs = 2, metavar = 'float float',
help = 'data range (min max) [auto]')
parser.add_option('--gap', '--transparent',
dest = 'gap',
type = 'float', metavar = 'float',
help = 'value to treat as transparent [%default]')
parser.add_option('-d','--dimension',
dest = 'dimension',
type = 'int', nargs = 3, metavar = ' '.join(['int']*3),
help = 'data dimension (x/y/z)')
parser.add_option('-s','--size',
dest = 'size',
type = 'float', nargs = 3, metavar = ' '.join(['float']*3),
help = 'box size (x/y/z)')
parser.add_option('-f','--defgrad',
dest = 'defgrad', metavar = 'string',
help = 'column label of deformation gradient [%default]')
parser.add_option('--scaling',
dest = 'scaling',
type = 'float', nargs = 3, metavar = ' '.join(['float']*3),
help = 'x/y/z scaling of displacement fluctuation [%default]')
parser.add_option('-z','--layer',
dest = 'z',
type = 'int', metavar = 'int',
help = 'index of z plane to plot [%default]')
parser.add_option('--color',
dest = 'color',
type = 'string', metavar = 'string',
help = 'color scheme')
parser.add_option('--invert',
dest = 'invert',
action = 'store_true',
help = 'invert color scheme')
parser.add_option('--abs',
dest = 'abs',
action = 'store_true',
help = 'magnitude of values')
parser.add_option('--log',
dest = 'log',
action = 'store_true',
help = 'log10 of values')
parser.add_option('-N','--pixelsize',
dest = 'pixelsize',
type = 'int', metavar = 'int',
help = 'pixels per cell edge')
parser.add_option('--show',
dest = 'show',
action = 'store_true',
help = 'show resulting image')
parser.set_defaults(label = None,
range = [0.0,0.0],
dimension = [],
size = [],
z = 1,
abs = False,
log = False,
defgrad = 'f',
scaling = [1.,1.,1.],
color = "gray",
invert = False,
pixelsize = 1,
show = False,
)
(options,filenames) = parser.parse_args()
if filenames == []: filenames = [None]
options.size = np.array(options.size)
options.dimension = np.array(options.dimension)
options.range = np.array(options.range)
if options.z > 0: options.z -= 1 # adjust to 0-based indexing
# --- color palette ---------------------------------------------------------------------------------
theMap = damask.Colormap(predefined=options.color)
if options.invert: theMap = theMap.invert()
theColors = np.uint8(np.array(theMap.export(format='list',steps=256))*255)
# --- loop over input files -------------------------------------------------------------------------
for name in filenames:
try:
table = damask.ASCIItable(name = name, labeled = options.label is not None, readonly = True)
except IOError:
continue
damask.util.report(scriptName,name)
# ------------------------------------------ read header ------------------------------------------
table.head_read()
# --------------- figure out columns to process ---------------------------------------------------
errors = []
if table.label_dimension(options.label) != 1:
errors.append('no scalar data ({}) found.'.format(options.label))
if table.label_dimension(options.defgrad) != 9:
errors.append('no deformation gradient tensor (1..9_{}) found.'.format(options.defgrad))
if errors != []:
table.croak('\n'.join(errors)+'\n')
table.close(dismiss = True)
continue
table.data_readArray([options.label,options.defgrad])
data = table.data[:,0 ].transpose().reshape( list(options.dimension),order='F')
F = table.data[:,1:10].transpose().reshape([3,3]+list(options.dimension),order='F')
if options.abs: data = np.abs(data)
if options.log: data = np.log10(data)
if np.all(options.range == 0.0): options.range = np.array([data.min(),data.max()])
elif options.log: options.range = np.log10(options.range)
data = ( data - options.range.min()) / \
(options.range.max() - options.range.min()) # data scaled to fraction of range
data = np.clip(data,0.0,1.0) # cut off outliers (should be none)
# ---------------- calculate coordinates -----------------------------------------------------------
Favg = damask.core.math.tensorAvg(F)
centroids = damask.core.mesh.deformedCoordsFFT(options.size,F,Favg,options.scaling)
nodes = damask.core.mesh.nodesAroundCentres(options.size,Favg,centroids)
boundingBox = np.array([ \
[np.amin(nodes[0,:,:,options.z]),np.amin(nodes[1,:,:,options.z]),np.amin(nodes[2,:,:,options.z])],
[np.amax(nodes[0,:,:,options.z]),np.amax(nodes[1,:,:,options.z]),np.amax(nodes[2,:,:,options.z])],
]) # find x-y bounding box for given z layer
nodes -= boundingBox[0].repeat(np.prod(options.dimension+1)).reshape([3]+list(options.dimension+1))
nodes *= (options.pixelsize*options.dimension/options.size).repeat(np.prod(options.dimension+1)).\
reshape([3]+list(options.dimension+1))
imagesize = (options.pixelsize*(boundingBox[1]-boundingBox[0])* # determine image size from number of
options.dimension/options.size)[:2].astype('i') # cells in overall bounding box
im = Image.new('RGBA',imagesize)
draw = ImageDraw.Draw(im)
for y in range(options.dimension[1]):
for x in range(options.dimension[0]):
draw.polygon([nodes[0,x ,y ,options.z],
nodes[1,x ,y ,options.z],
nodes[0,x+1,y ,options.z],
nodes[1,x+1,y ,options.z],
nodes[0,x+1,y+1,options.z],
nodes[1,x+1,y+1,options.z],
nodes[0,x ,y+1,options.z],
nodes[1,x ,y+1,options.z],
],
fill = tuple(theColors[int(255*data[x,y,options.z])],
0 if data[x,y,options.z] == options.gap else 255),
outline = None)
# ------------------------------------------ output result -----------------------------------------
im.save(os.path.splitext(name)[0]+ \
('_'+options.label if options.label else '')+ \
'.png' if name else sys.stdout,
format = "PNG")
table.close() # close ASCII table
if options.show: im.show()

0
processing/pre/seeds_fromDistribution.py → processing/legacy/seeds_fromDistribution.py
View File

454
processing/legacy/vtk2ang.py
View File

@ -1,454 +0,0 @@
#!/usr/bin/env python3
import os
import sys
from optparse import OptionParser
import string
import math
import numpy as np
import vtk
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
# -----------------------------
def getHeader(filename,sizeFastIndex,sizeSlowIndex,stepsize):
"""Returns header for ang file step size in micrometer."""
return '\n'.join([
'# TEM_PIXperUM 1.000000',
'# x-star 1.000000',
'# y-star 1.000000',
'# z-star 1.000000',
'# WorkingDistance 18.000000',
'#',
'# Phase 1',
'# MaterialName XX',
'# Formula XX',
'# Info',
'# Symmetry 43',
'# LatticeConstants 2.870 2.870 2.870 90.000 90.000 90.000',
'# NumberFamilies 1',
'# hklFamilies 1 1 0 1 0.000000 1',
'# Categories 0 0 0 0 0 ',
'#',
'# GRID: SqrGrid',
'# XSTEP: ' + str(stepsize*1e6),
'# YSTEP: ' + str(stepsize*1e6),
'# NCOLS_ODD: ' + str(sizeFastIndex),
'# NCOLS_EVEN: ' + str(sizeFastIndex),
'# NROWS: ' + str(sizeSlowIndex),
'#',
'# OPERATOR: ' + string.replace('$Id$','\n','\\n'),
'#',
'# SAMPLEID: {}'.format(filename),
'#',
'# SCANID: ',
'#',
]) + '\n'
# -----------------------------
def positiveRadians(angle):
"""Returns positive angle in radians from angle in degrees."""
return math.radians(float(angle)) % (2.*math.pi)
# -----------------------------
def getDataLine(angles,x,y,validData=True):
"""
Returns string of one line in ang file.
Convention in ang file: y coordinate comes first and is fastest index
positions in micrometer.
"""
info = {True: (9999.9, 1.0, 0,99999,0.0),
False: ( -1.0,-1.0,-1, -1,1.0)}
return '%9.5f %9.5f %9.5f %12.5f %12.5f %6.1f %6.3f %2i %6i %6.3f \n'\
%(tuple(map(positiveRadians,angles))+(y*1e6,x*1e6)+info[validData])
# --------------------------------------------------------------------
# MAIN FUNCTION STARTS HERE
# --------------------------------------------------------------------
parser = OptionParser(usage='%prog options [file[s]]', description = """
Builds a ang files from a vtk file.
""", version = scriptID)
parser.add_option('--disp','--displacement',dest='dispLabel',
metavar ='string',
help='label of displacements [%default]')
parser.add_option('--euler', dest='eulerLabel', nargs=3,
metavar ='string string string',
help='labels of euler angles [%default]')
parser.add_option('-n','--normal', dest='normal', type='float', nargs=3,
metavar ='float float float',
help='normal of slices in direction of increasing slice numbers [%default]')
parser.add_option('-u','--up', dest='up', type='float', nargs=3,
metavar ='float float float',
help='up direction of slices [%default]')
parser.add_option('-i','--slices', dest='Nslices', type='int',
metavar ='int',
help='number of slices [%default]')
parser.add_option('-d','--distance', dest='distance', type='float',
metavar ='float',
help='slice distance [%default]')
parser.add_option('-s','--scale', dest='scale', type='float',
metavar ='float',
help='scale length from vtk file [%default]')
parser.add_option('-r','--resolution', dest='resolution', type='float',
metavar ='float',
help='scaling factor for resolution [%default]')
parser.add_option('--verbose', dest='verbose', action='store_true',
help='verbose mode')
parser.add_option('--visualize', dest='visualize', action='store_true',
help='visualize geometry')
parser.set_defaults(dispLabel = 'displacement')
parser.set_defaults(eulerLabel = ['1_1_eulerangles','1_2_eulerangles','1_3_eulerangles'])
parser.set_defaults(normal = [0.0,0.0,-1.0])
parser.set_defaults(up = [0.0,1.0,0.0])
parser.set_defaults(Nslices = 1)
parser.set_defaults(distance = 0.0)
parser.set_defaults(scale = 1.0)
parser.set_defaults(resolution = 1.0)
parser.set_defaults(dispScaling = 1.0)
parser.set_defaults(verbose = False)
parser.set_defaults(visualize = False)
(options,filenames) = parser.parse_args()
#--- SANITY CHECKS
# check for valid filenames
for filename in filenames:
if not os.path.exists(filename):
parser.error('file "%s" does not exist'%filename)
if not os.path.splitext(filename)[1] == '.vtk':
parser.error('"%s": need vtk file'%filename)
# check for othogonality of normal and up vector
if np.dot(np.array(options.normal),np.array(options.up)) > 1e-3:
parser.error('normal vector and up vector have to be orthogonal')
#--- ITERATE OVER FILES AND PROCESS THEM
for filename in filenames:
if options.verbose: sys.stdout.write("\nREADING VTK FILE\n")
# Read the source file
reader = vtk.vtkUnstructuredGridReader()
reader.SetFileName(filename)
reader.ReadAllScalarsOn()
reader.ReadAllVectorsOn()
reader.Update()
undeformedMesh = reader.GetOutput()
# Get euler angles from cell data
if options.verbose: sys.stdout.write("\nGETTING EULER ANGLES\n")
angles = {}
for i in range(undeformedMesh.GetPointData().GetNumberOfArrays()):
scalarName = undeformedMesh.GetPointData().GetArrayName(i)
if scalarName in options.eulerLabel:
angles[scalarName] = undeformedMesh.GetCellData().GetScalars(scalarName)
if options.verbose: sys.stdout.write(" found scalar with name %s\n"%scalarName)
if len(angles) < 3: # found data for all three euler angles?
for label in options.eulerLabel:
if label not in angles.keys():
parser.error('Could not find scalar data with name %s'%label)
# Get deformed mesh
if options.verbose: sys.stdout.write("\nDEFORM MESH\n")
warpVector = vtk.vtkWarpVector()
undeformedMesh.GetPointData().SetActiveVectors(options.dispLabel)
warpVector.SetInputData(undeformedMesh)
warpVector.Update()
deformedMesh = warpVector.GetOutput()
box = deformedMesh.GetBounds() # bounding box in mesh system
if options.verbose:
sys.stdout.write(" bounding box in lab system\n")
sys.stdout.write(" x (% .8f % .8f)\n"%(box[0],box[1]))
sys.stdout.write(" y (% .8f % .8f)\n"%(box[2],box[3]))
sys.stdout.write(" z (% .8f % .8f)\n"%(box[4],box[5]))
# Get cell centers of deformed mesh (position of ips)
if options.verbose: sys.stdout.write("\nGETTING CELL CENTERS OF DEFORMED MESH\n")
cellCenter = vtk.vtkCellCenters()
cellCenter.SetVertexCells(0) # do not generate vertex cells, just points
cellCenter.SetInputData(deformedMesh)
cellCenter.Update()
meshIPs = cellCenter.GetOutput()
# Get outer surface of deformed mesh
if options.verbose: sys.stdout.write("\nGETTING OUTER SURFACE OF DEFORMED MESH\n")
surfaceFilter = vtk.vtkDataSetSurfaceFilter()
surfaceFilter.SetInputData(deformedMesh)
surfaceFilter.Update()
surface = surfaceFilter.GetOutput()
# Get coordinate system for ang files
# z-vector is normal to slices
# x-vector corresponds to the up-direction
# "R" rotates coordinates from the mesh system into the TSL system
if options.verbose: sys.stdout.write("\nGETTING COORDINATE SYSTEM FOR ANG FILES\n")
z = np.array(options.normal,dtype='float')
z = z / np.linalg.norm(z)
x = np.array(options.up,dtype='float')
x = x / np.linalg.norm(x)
y = np.cross(z,x)
R = np.array([x,y,z])
if options.verbose:
sys.stdout.write(" axis (x: up direction, z: slice normal)\n")
sys.stdout.write(" x (% .8f % .8f % .8f)\n"%tuple(x))
sys.stdout.write(" y (% .8f % .8f % .8f)\n"%tuple(y))
sys.stdout.write(" z (% .8f % .8f % .8f)\n"%tuple(z))
# Get bounding box in rotated system (x,y,z)
if options.verbose: sys.stdout.write("\nGETTING BOUNDING BOX IN ROTATED SYSTEM\n")
rotatedbox = [[np.inf,-np.inf] for i in range(3)] # bounding box in rotated TSL system
for n in range(8): # loop over eight vertices of mesh bounding box
vert = np.array([box[0+(n//1)%2],
box[2+(n//2)%2],
box[4+(n//4)%2]]) # vertex in mesh system
rotatedvert = np.dot(R,vert) # vertex in rotated system
for i in range(3):
rotatedbox[i][0] = min(rotatedbox[i][0],rotatedvert[i])
rotatedbox[i][1] = max(rotatedbox[i][1],rotatedvert[i])
if options.verbose:
sys.stdout.write(" bounding box in rotated system\n")
sys.stdout.write(" x (% .8f % .8f)\n"%tuple(rotatedbox[0]))
sys.stdout.write(" y (% .8f % .8f)\n"%tuple(rotatedbox[1]))
sys.stdout.write(" z (% .8f % .8f)\n"%tuple(rotatedbox[2]))
# Correct bounding box so that a multiplicity of the resolution fits into it
# and get number of points and extent in each (rotated) axis direction
if options.verbose: sys.stdout.write("\nCORRECTING EXTENT OF BOUNDING BOX IN ROTATED SYSTEM\n")
correction = []
Npoints = []
extent = [rotatedbox[i][1] - rotatedbox[i][0] for i in range(3)]
for i in range(2):
Npoints.extend([int(math.ceil(extent[i] / options.resolution))])
correction.extend([float(Npoints[i]) * options.resolution - extent[i]])
if options.distance > 0.0:
Npoints.extend([int(math.ceil(extent[2] / options.distance))])
correction.extend([float(Npoints[2]) * options.distance - extent[2]])
else:
Npoints.extend([options.Nslices])
correction.extend([0.0])
options.distance = extent[2] / float(options.Nslices)
for i in range(3):
rotatedbox[i][0] -= 0.5 * correction[i]
rotatedbox[i][1] += 0.5 * correction[i]
extent[i] = rotatedbox[i][1] - rotatedbox[i][0]
NpointsPerSlice = Npoints[0] * Npoints[1]
totalNpoints = NpointsPerSlice * Npoints[2]
if options.verbose:
sys.stdout.write(" corrected bounding box in rotated system\n")
sys.stdout.write(" x (% .8f % .8f)\n"%tuple(rotatedbox[0]))
sys.stdout.write(" y (% .8f % .8f)\n"%tuple(rotatedbox[1]))
sys.stdout.write(" z (% .8f % .8f)\n"%tuple(rotatedbox[2]))
# Generate new regular point grid for ang files
# Use "polydata" object with points as single vertices
# beware of TSL convention: y direction is fastest index
if options.verbose: sys.stdout.write("\nGENERATING POINTS FOR POINT GRID")
points = vtk.vtkPoints()
for k in range(Npoints[2]):
for j in range(Npoints[0]):
for i in range(Npoints[1]): # y is fastest index
rotatedpoint = np.array([rotatedbox[0][0] + (float(j) + 0.5) * options.resolution,
rotatedbox[1][0] + (float(i) + 0.5) * options.resolution,
rotatedbox[2][0] + (float(k) + 0.5) * options.distance ]) # point in rotated system
point = np.dot(R.T,rotatedpoint) # point in mesh system
points.InsertNextPoint(list(point))
if options.verbose:
sys.stdout.write("\rGENERATING POINTS FOR POINT GRID %d%%" %(100*(Npoints[1]*(k*Npoints[0]+j)+i+1)/totalNpoints))
sys.stdout.flush()
if options.verbose:
sys.stdout.write("\n number of slices: %i\n"%Npoints[2])
sys.stdout.write(" slice spacing: %.8f\n"%options.distance)
if Npoints[2] > 1:
sys.stdout.write(" number of points per slice: %i = %i rows * %i points in row\n"%(NpointsPerSlice,Npoints[0],Npoints[1]))
sys.stdout.write(" grid resolution: %.8f\n"%options.resolution)
if options.verbose: sys.stdout.write("\nGENERATING VERTICES FOR POINT GRID")
vertices = vtk.vtkCellArray()
for i in range(totalNpoints):
vertex = vtk.vtkVertex()
vertex.GetPointIds().SetId(0,i) # each vertex consists of exactly one (index 0) point with ID "i"
vertices.InsertNextCell(vertex)
if options.verbose:
sys.stdout.write("\rGENERATING VERTICES FOR POINT GRID %d%%" %(100*(i+1)/totalNpoints))
sys.stdout.flush()
if options.verbose: sys.stdout.write("\n\nGENERATING POINT GRID\n")
pointgrid = vtk.vtkPolyData()
pointgrid.SetPoints(points)
pointgrid.SetVerts(vertices)
pointgrid.Update()
# Find out which points reside inside mesh geometry
if options.verbose: sys.stdout.write("\nIDENTIFYING POINTS INSIDE MESH GEOMETRY\n")
enclosedPoints = vtk.vtkSelectEnclosedPoints()
enclosedPoints.SetSurface(surface)
enclosedPoints.SetInput(pointgrid)
enclosedPoints.Update()
# Build kdtree from mesh IPs and match mesh IPs to point grid
if options.verbose: sys.stdout.write("\nBUILDING MAPPING OF GRID POINTS")
kdTree = vtk.vtkKdTree()
kdTree.BuildLocatorFromPoints(meshIPs.GetPoints())
gridToMesh = []
ids = vtk.vtkIdList()
NenclosedPoints = 0
for i in range(pointgrid.GetNumberOfPoints()):
gridToMesh.append([])
if enclosedPoints.IsInside(i):
NenclosedPoints += 1
# here one could use faster(?) "FindClosestPoint" if only first nearest neighbor required
kdTree.FindClosestNPoints(1,pointgrid.GetPoint(i),ids)
for j in range(ids.GetNumberOfIds()):
gridToMesh[-1].extend([ids.GetId(j)])
if options.verbose:
sys.stdout.write("\rBUILDING MAPPING OF GRID POINTS %d%%" %(100*(i+1)/totalNpoints))
sys.stdout.flush()
if options.verbose:
sys.stdout.write("\n Number of points inside mesh geometry %i\n"%NenclosedPoints)
sys.stdout.write(" Number of points outside mesh geometry %i\n"%(totalNpoints - NenclosedPoints))
# ITERATE OVER SLICES AND CREATE ANG FILE
if options.verbose:
sys.stdout.write("\nWRITING OUT ANG FILES\n")
sys.stdout.write(" scaling all length with %f\n"%options.scale)
x0,y0,z0 = np.dot(R,pointgrid.GetPoint(0)) # first point on slice defines origin
for sliceN in range(Npoints[2]):
# Open file and write header
angfilename = eval('"'+eval("'%%s_slice%%0%ii.ang'%(math.log10(Npoints[2])+1)")+'"%(os.path.splitext(filename)[0],sliceN+1)')
with open(angfilename,'w') as angfile:
if options.verbose: sys.stdout.write(" %s\n"%angfilename)
angfile.write(getHeader(filename,Npoints[1],Npoints[0],options.resolution*options.scale))
for i in range(sliceN*NpointsPerSlice,(sliceN+1)*NpointsPerSlice): # Iterate over points on slice
# Get euler angles of closest IDs
if enclosedPoints.IsInside(i):
phi = []
for j in range(len(gridToMesh[i])):
IP = gridToMesh[i][j]
phi.append([])
for k in range(3):
phi[-1].extend([angles[options.eulerLabel[k]].GetValue(IP)])
else:
phi = [[720,720,720]] # fake angles
# Interpolate Euler angle
# NOT YET IMPLEMENTED, simply take the nearest neighbors values
interpolatedPhi = phi[0]
# write data to ang file
x,y,z = np.dot(R,pointgrid.GetPoint(i)) # point in rotated TSL system
x -= x0 # first point on slice defines origin
y -= y0 # first point on slice defines origin
x *= options.scale
y *= options.scale
angfile.write(getDataLine(interpolatedPhi,x,y,enclosedPoints.IsInside(i)))
# Visualize slices
if options.visualize:
meshMapper = vtk.vtkDataSetMapper()
meshMapper.SetInput(surface)
meshMapper.ScalarVisibilityOff() # do not use scalar data for coloring
meshActor = vtk.vtkActor()
meshActor.SetMapper(meshMapper)
meshActor.GetProperty().SetOpacity(0.2)
meshActor.GetProperty().SetColor(1.0,1.0,0)
meshActor.GetProperty().BackfaceCullingOn()
boxpoints = vtk.vtkPoints()
for n in range(8):
P = [rotatedbox[0][(n/1)%2],
rotatedbox[1][(n/2)%2],
rotatedbox[2][(n/4)%2]]
boxpoints.InsertNextPoint(list(np.dot(R.T,np.array(P))))
box = vtk.vtkHexahedron()
for n,i in enumerate([0,1,3,2,4,5,7,6]):
box.GetPointIds().SetId(n,i)
boxgrid = vtk.vtkUnstructuredGrid()
boxgrid.SetPoints(boxpoints)
boxgrid.InsertNextCell(box.GetCellType(), box.GetPointIds())
boxsurfaceFilter = vtk.vtkDataSetSurfaceFilter()
boxsurfaceFilter.SetInput(boxgrid)
boxsurfaceFilter.Update()
boxsurface = boxsurfaceFilter.GetOutput()
boxMapper = vtk.vtkDataSetMapper()
boxMapper.SetInput(boxsurface)
boxActor = vtk.vtkActor()
boxActor.SetMapper(boxMapper)
boxActor.GetProperty().SetLineWidth(2.0)
boxActor.GetProperty().SetRepresentationToWireframe()
gridMapper = vtk.vtkDataSetMapper()
gridMapper.SetInput(pointgrid)
gridActor = vtk.vtkActor()
gridActor.SetMapper(gridMapper)
gridActor.GetProperty().SetColor(0,0,0)
gridActor.GetProperty().SetPointSize(3)
renderer = vtk.vtkRenderer()
renderWindow = vtk.vtkRenderWindow()
renderWindow.FullScreenOn()
renderWindow.AddRenderer(renderer)
renderWindowInteractor = vtk.vtkRenderWindowInteractor()
renderWindowInteractor.SetRenderWindow(renderWindow)
renderer.AddActor(meshActor)
renderer.AddActor(boxActor)
renderer.AddActor(gridActor)
renderer.SetBackground(1,1,1)
renderWindow.Render()
renderWindowInteractor.SetInteractorStyle(vtk.vtkInteractorStyleTrackballCamera())
renderWindowInteractor.Start()

63
processing/post/DADF5_postResults.py
View File

@ -1,63 +0,0 @@
#!/usr/bin/env python3
import os
import argparse
import numpy as np
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
parser = argparse.ArgumentParser()
parser.add_argument('filenames', nargs='+',
help='DADF5 files')
parser.add_argument('-d','--dir', dest='dir',default='postProc',metavar='string',
help='name of subdirectory relative to the location of the DADF5 file to hold output')
parser.add_argument('--mat', nargs='+',
help='labels for homogenization',dest='mat')
parser.add_argument('--con', nargs='+',
help='labels for phase',dest='con')
options = parser.parse_args()
if options.mat is None: options.mat=[]
if options.con is None: options.con=[]
for filename in options.filenames:
results = damask.Result(filename)
if not results.structured: continue
coords = damask.grid_filters.coordinates0_point(results.cells,results.size,results.origin).reshape(-1,3,order='F')
N_digits = int(np.floor(np.log10(int(results.increments[-1][10:]))))+1
N_digits = 5 # hack to keep test intact
for inc in damask.util.show_progress(results.iterate('increments'),len(results.increments)):
table = damask.Table(np.ones(np.product(results.cells),dtype=int)*int(inc[10:]),{'inc':(1,)})\
.add('pos',coords.reshape(-1,3))
results.view('homogenizations',False)
results.view('phases',True)
for label in options.con:
x = results.get_dataset_location(label)
if len(x) != 0:
table = table.add(label,results.read_dataset(x,0,plain=True).reshape(results.cells.prod(),-1))
results.view('phases',False)
results.view('homogenizations',True)
for label in options.mat:
x = results.get_dataset_location(label)
if len(x) != 0:
table = table.add(label,results.read_dataset(x,0,plain=True).reshape(results.cells.prod(),-1))
dirname = os.path.abspath(os.path.join(os.path.dirname(filename),options.dir))
if not os.path.isdir(dirname):
os.mkdir(dirname,0o755)
file_out = '{}_inc{}.txt'.format(os.path.splitext(os.path.split(filename)[-1])[0],
inc[10:].zfill(N_digits))
table.save(os.path.join(dirname,file_out),legacy=True)

71
processing/pre/geom_fromDREAM3D.py
View File

@ -1,71 +0,0 @@
#!/usr/bin/env python3
import os
from optparse import OptionParser
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
#--------------------------------------------------------------------------------------------------
# MAIN
#--------------------------------------------------------------------------------------------------
parser = OptionParser(usage='%prog options [DREAM.3Dfile(s)]', description = """
Converts DREAM.3D file. Input can be cell data (direct pointwise takeover) or grain data (individual
grains are segmented). Requires orientation data as quaternion.
""", version = scriptID)
parser.add_option('-b','--basegroup',
dest = 'basegroup',
metavar = 'string',
help = 'name of the group in "DataContainers" containing the pointwise (and, if applicable grain average) data')
parser.add_option('-p','--pointwise',
dest = 'pointwise',
metavar = 'string',
help = 'name of the group in "DataContainers/<basegroup>" containing pointwise data [%default]')
parser.add_option('-a','--average',
dest = 'average',
metavar = 'string',
help = 'name of the group in "DataContainers</basegroup>" containing grain average data. '\
+ 'Leave empty for pointwise data')
parser.add_option('--phase',
dest = 'phase',
type = 'string',
metavar = 'string',
help = 'name of the dataset containing pointwise/average phase IDs [%default]')
parser.add_option('--microstructure',
dest = 'microstructure',
type = 'string',
metavar = 'string',
help = 'name of the dataset connecting pointwise and average data [%default]')
parser.add_option('-q', '--quaternion',
dest = 'quaternion',
type = 'string',
metavar='string',
help = 'name of the dataset containing pointwise/average orientation as quaternion [%default]')
parser.set_defaults(pointwise = 'CellData',
quaternion = 'Quats',
phase = 'Phases',
microstructure = 'FeatureIds',
)
(options, filenames) = parser.parse_args()
if options.basegroup is None:
parser.error('No base group selected')
if filenames == []: parser.error('no input file specified.')
for name in filenames:
damask.util.report(scriptName,name)
geom = damask.Grid.load_DREAM3D(name,'FeatureIds')
damask.util.croak(geom)
geom.save_ASCII(os.path.splitext(name)[0]+'.geom')

141
processing/pre/geom_fromOsteonGeometry.py
View File

@ -1,141 +0,0 @@
#!/usr/bin/env python3
import os
import sys
from optparse import OptionParser
import numpy as np
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
parser = OptionParser(usage='%prog options [geomfile]', description = """
Generate description of an osteon enclosing the Harvesian canal and separated by interstitial tissue.
The osteon phase is lamellar with a twisted plywood structure.
Its fiber orientation is oscillating by +/- amplitude within one period.
""", version = scriptID)
parser.add_option('-g', '--grid',
dest='grid', type='int',
nargs=2, metavar = 'int int',
help='a,b grid of hexahedral box [%default]')
parser.add_option('-s', '--size',
dest='size',
type='float', nargs=2, metavar = 'float float',
help='x,y size of hexahedral box [%default]')
parser.add_option('-c', '--canal',
dest='canal',
type='float', metavar = 'float',
help='Haversian canal radius [%default]')
parser.add_option('-o', '--osteon',
dest='osteon',
type='float', metavar = 'float',
help='horizontal osteon radius [%default]')
parser.add_option('-l', '--lamella',
dest='period',
type='float', metavar = 'float',
help='lamella width [%default]')
parser.add_option('-a', '--amplitude',
dest='amplitude',
type='float', metavar = 'float',
help='amplitude of twisted plywood wiggle in deg [%default]')
parser.add_option( '--aspect',
dest='aspect',
type='float', metavar = 'float',
help='vertical/horizontal osteon aspect ratio [%default]')
parser.add_option('-w', '--omega',
dest='omega',
type='float', metavar = 'float',
help='rotation angle around normal of osteon [%default]')
parser.set_defaults(canal = 25e-6,
osteon = 100e-6,
aspect = 1.0,
omega = 0.0,
period = 5e-6,
amplitude = 60,
size = (300e-6,300e-6),
grid = (512,512),
)
(options,filename) = parser.parse_args()
name = None if filename == [] else filename[0]
damask.util.report(scriptName,name)
omega = np.deg2rad(options.omega)
rotation = np.array([[ np.cos(omega),np.sin(omega),],
[-np.sin(omega),np.cos(omega),]])
grid = np.array(options.grid,'i')
size = np.array(options.size,'d')
X0,Y0 = np.meshgrid(size[0]/grid[0] * (np.arange(grid[0]) - grid[0]/2 + 0.5),
size[1]/grid[0] * (np.arange(grid[1]) - grid[1]/2 + 0.5), indexing='ij')
X = X0*rotation[0,0] + Y0*rotation[0,1] # rotate by omega
Y = X0*rotation[1,0] + Y0*rotation[1,1] # rotate by omega
radius = np.sqrt(X*X + Y*Y/options.aspect**2.0)
alpha = np.degrees(np.arctan2(Y/options.aspect,X))
beta = options.amplitude*np.sin(2.0*np.pi*(radius-options.canal)/options.period)
microstructure = np.where(radius < float(options.canal), 1,0) \
+ np.where(radius > float(options.osteon),2,0)
# extend to 3D
size = np.append(size,np.min(size/grid))
grid = np.append(grid,1)
microstructure = microstructure.reshape(microstructure.shape+(1,))
Alpha = np.zeros(grid[0]*grid[1],'d')
Beta = np.zeros(grid[0]*grid[1],'d')
i = 3
for y in range(grid[1]):
for x in range(grid[0]):
if microstructure[x,y] == 0:
microstructure[x,y] = i
Alpha[i] = alpha[x,y]
Beta [i] = beta [x,y]
i+=1
config_header = ['<texture>',
'[canal]',
'[interstitial]'
]
for i in range(3,np.max(microstructure)):
config_header += ['[Point{}]'.format(i-2),
'(gauss)\tphi1 {:.2f}\tPhi {:.2f}\tphi2 0'.format(Alpha[i],Beta[i])
]
config_header = ['<microstructure>',
'[canal]',
'(constituent)\tphase 1\ttexture 1\tfraction 1.0',
'[interstitial]',
'(constituent)\tphase 2\ttexture 2\tfraction 1.0'
]
for i in range(3,np.max(microstructure)):
config_header += ['[Point{}]'.format(i-2),
'(constituent)\tphase 3\ttexture {}\tfraction 1.0'.format(i)
]
header = [scriptID + ' ' + ' '.join(sys.argv[1:])]\
+ config_header
geom = damask.Grid(microstructure.reshape(grid),
size,-size/2,
comments=header)
damask.util.croak(geom)
geom.save_ASCII(sys.stdout if name is None else name)

96
processing/pre/seeds_fromPokes.py
View File

@ -1,96 +0,0 @@
#!/usr/bin/env python3
import os
import sys
from io import StringIO
from optparse import OptionParser
import numpy as np
import damask
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version])
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
parser = OptionParser(usage='%prog options [file[s]]', description = """
Create seeds file by poking at 45 degree through given geom file.
Mimics APS Beamline 34-ID-E DAXM poking.
""", version = scriptID)
parser.add_option('-N', '--points',
dest = 'N',
type = 'int', metavar = 'int',
help = 'number of poking locations [%default]')
parser.add_option('-b', '--box',
dest = 'box',
type = 'float', nargs = 6, metavar = ' '.join(['float']*6),
help = 'bounding box as fraction in x, y, and z directions')
parser.add_option('-x',
action = 'store_true',
dest = 'x',
help = 'poke 45 deg along x')
parser.add_option('-y',
action = 'store_true',
dest = 'y',
help = 'poke 45 deg along y')
parser.set_defaults(x = False,
y = False,
box = [0.0,1.0,0.0,1.0,0.0,1.0],
N = 16,
)
(options,filenames) = parser.parse_args()
if filenames == []: filenames = [None]
options.box = np.array(options.box).reshape(3,2)
for name in filenames:
damask.util.report(scriptName,name)
geom = damask.Grid.load_ASCII(StringIO(''.join(sys.stdin.read())) if name is None else name)
offset =(np.amin(options.box, axis=1)*geom.cells/geom.size).astype(int)
box = np.amax(options.box, axis=1) \
- np.amin(options.box, axis=1)
Nx = int(options.N/np.sqrt(options.N*geom.size[1]*box[1]/geom.size[0]/box[0]))
Ny = int(options.N/np.sqrt(options.N*geom.size[0]*box[0]/geom.size[1]/box[1]))
Nz = int(box[2]*geom.cells[2])
damask.util.croak('poking {} x {} x {} in box {} {} {}...'.format(Nx,Ny,Nz,*box))
seeds = np.zeros((Nx*Ny*Nz,4))
g = np.zeros(3,dtype=np.int)
n = 0
for i in range(Nx):
for j in range(Ny):
g[0] = round((i+0.5)*box[0]*geom.cells[0]/Nx-0.5)+offset[0]
g[1] = round((j+0.5)*box[1]*geom.cells[1]/Ny-0.5)+offset[1]
for k in range(Nz):
g[2] = k + offset[2]
g %= geom.cells
seeds[n,0:3] = (g+0.5)/geom.cells # normalize coordinates to box
seeds[n, 3] = geom.material[g[0],g[1],g[2]]
if options.x: g[0] += 1
if options.y: g[1] += 1
n += 1
comments = geom.comments \
+ [scriptID + ' ' + ' '.join(sys.argv[1:]),
'poking\ta {}\tb {}\tc {}'.format(Nx,Ny,Nz),
'grid\ta {}\tb {}\tc {}'.format(*geom.cells),
'size\tx {}\ty {}\tz {}'.format(*geom.size),
'origin\tx {}\ty {}\tz {}'.format(*geom.origin),
]
table = damask.Table(seeds,{'pos':(3,),'material':(1,)},comments)
table.set('material',table.get('material').astype(np.int))\
.save(sys.stdout if name is None else \
os.path.splitext(name)[0]+f'_poked_{options.N}.seeds',legacy=True)