switched to new keywords "grid" and "size", still reads the old "resolution" and "dimension" but writes the new ones

This commit is contained in:
Martin Diehl 2013-04-12 10:27:05 +00:00
parent 964b0d7839
commit 0b1033d63a
4 changed files with 185 additions and 160 deletions

View File

@ -1,5 +1,6 @@
#!/usr/bin/env python #!/usr/bin/env python
# -*- coding: utf-8 -*- # -*- coding: utf-8 -*-
import os,sys,math,string,re,numpy, damask import os,sys,math,string,re,numpy, damask
from optparse import OptionParser, OptionGroup, Option, SUPPRESS_HELP from optparse import OptionParser, OptionGroup, Option, SUPPRESS_HELP
@ -23,16 +24,14 @@ class extendedOption(Option):
Option.take_action(self, action, dest, opt, value, values, parser) Option.take_action(self, action, dest, opt, value, values, parser)
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# MAIN # MAIN
# -------------------------------------------------------------------- # --------------------------------------------------------------------
identifiers = { identifiers = {
'resolution': ['a','b','c'], 'grid': ['a','b','c'],
} }
mappings = { mappings = {
'resolution': lambda x: int(x), 'grid': lambda x: int(x),
'grains': lambda x: int(x), 'grains': lambda x: int(x),
} }
@ -41,24 +40,24 @@ Generate geometry description and material configuration by standard Voronoi tes
""" + string.replace('$Id$','\n','\\n') """ + string.replace('$Id$','\n','\\n')
) )
parser.add_option('-r', '--resolution', dest='resolution', type='int', nargs = 3, \ parser.add_option('-g', '--grid', dest='grid', type='int', nargs = 3, \
help='a,b,c resolution of periodic box') help='a,b,c grid of hexahedral box [from seed file]')
parser.add_option('-d', '--dimension', dest='dimension', type='float', nargs = 3, \ parser.add_option('-s', '--size', dest='size', type='float', nargs = 3, \
help='x,y,z dimension of periodic box') help='x,y,z size of hexahedral box [1.0 along largest grid point number]')
parser.add_option('--homogenization', dest='homogenization', type='int', \ parser.add_option('--homogenization', dest='homogenization', type='int', \
help='homogenization index to be used') help='homogenization index to be used [%default]')
parser.add_option('--phase', dest='phase', type='int', \ parser.add_option('--phase', dest='phase', type='int', \
help='phase index to be used') help='phase index to be used [%default]')
parser.add_option('--crystallite', dest='crystallite', type='int', \ parser.add_option('--crystallite', dest='crystallite', type='int', \
help='crystallite index to be used') help='crystallite index to be used [%default]')
parser.add_option('-c', '--configuration', dest='config', action='store_true', \ parser.add_option('-c', '--configuration', dest='config', action='store_true', \
help='output material configuration') help='output material configuration [%default]')
parser.add_option('-2', '--twodimensional', dest='twoD', action='store_true', \ parser.add_option('-2', '--twodimensional', dest='twoD', action='store_true', \
help='output geom file with two-dimensional data arrangement') help='output geom file with two-dimensional data arrangement [%default]')
parser.set_defaults(resolution = [0,0,0]) parser.set_defaults(grid = [0,0,0])
parser.set_defaults(dimension = [0.0,0.0,0.0]) parser.set_defaults(size = [0.0,0.0,0.0])
parser.set_defaults(homogenization = 1) parser.set_defaults(homogenization = 1)
parser.set_defaults(phase = 1) parser.set_defaults(phase = 1)
parser.set_defaults(crystallite = 1) parser.set_defaults(crystallite = 1)
@ -105,8 +104,8 @@ for file in files:
file['input'].close() file['input'].close()
info = {'grains': 0, info = {'grains': 0,
'resolution': numpy.array([0,0,0]), 'grid': numpy.array([0,0,0]),
'dimension': numpy.array(options.dimension), 'size': numpy.array(options.size),
'origin': numpy.array([0.0,0.0,0.0]), 'origin': numpy.array([0.0,0.0,0.0]),
'homogenization': options.homogenization, 'homogenization': options.homogenization,
} }
@ -114,6 +113,7 @@ for file in files:
new_header = [] new_header = []
for header in headers: for header in headers:
headitems = map(str.lower,header.split()) headitems = map(str.lower,header.split())
if headitems[0] == 'resolution': headitems[0] = 'grid'
if headitems[0] in mappings.keys(): if headitems[0] in mappings.keys():
if headitems[0] in identifiers.keys(): if headitems[0] in identifiers.keys():
for i in xrange(len(identifiers[headitems[0]])): for i in xrange(len(identifiers[headitems[0]])):
@ -129,39 +129,26 @@ for file in files:
file['croak'].write('grain data not matching grain count...\n') file['croak'].write('grain data not matching grain count...\n')
info['grains'] = min(info['grains'],len(content)) info['grains'] = min(info['grains'],len(content))
if 0 not in options.resolution: # user-specified resolution if 0 not in options.grid: # user-specified grid
info['resolution'] = numpy.array(options.resolution) info['grid'] = numpy.array(options.grid)
if numpy.all(info['resolution'] == 0): if numpy.any(info['grid'] < 1):
file['croak'].write('no resolution info found.\n') file['croak'].write('no valid grid info found.\n')
continue continue
twoD = info['resolution'][2] < 2 twoD = info['grid'][2] < 2
for i in xrange(3): for i in xrange(3):
if info['dimension'][i] <= 0.0: # any invalid dimension? if info['size'][i] <= 0.0: # any invalid size?
info['dimension'][i] = float(info['resolution'][i])/max(info['resolution']) info['size'][i] = float(info['grid'][i])/max(info['grid'])
file['croak'].write('rescaling dimension %i...\n'%i) file['croak'].write('rescaling size %i...\n'%i)
file['croak'].write('grains: %i\n'%info['grains'] + \ file['croak'].write('grains to map: %i\n'%info['grains'] + \
'resolution: %s\n'%(' x '.join(map(str,info['resolution']))) + \ 'grid a b c: %s\n'%(' x '.join(map(str,info['grid']))) + \
'dimension: %s\n'%(' x '.join(map(str,info['dimension']))) + \ 'size x y z: %s\n'%(' x '.join(map(str,info['size']))) + \
'origin: %s\n'%(' : '.join(map(str,info['origin']))) + \ 'origin x y z: %s\n'%(' : '.join(map(str,info['origin']))) + \
'homogenization: %i\n'%info['homogenization']) 'homogenization: %i\n'%info['homogenization'])
new_header.append("resolution\ta %i\tb %i\tc %i\n"%(
info['resolution'][0],
info['resolution'][1],
info['resolution'][2],))
new_header.append("dimension\tx %f\ty %f\tz %f\n"%(
info['dimension'][0],
info['dimension'][1],
info['dimension'][2],))
new_header.append("origin\tx %f\ty %f\tz %f\n"%(
info['origin'][0],
info['origin'][1],
info['origin'][2],))
new_header.append("homogenization\t%i\n"%info['homogenization'])
# -------------------------------------- prepare data ---------------------------------- # -------------------------------------- prepare data ----------------------------------
@ -170,12 +157,12 @@ for file in files:
eulers = numpy.zeros((3,info['grains']),'d') eulers = numpy.zeros((3,info['grains']),'d')
for i in xrange(info['grains']): for i in xrange(info['grains']):
coords[:,i] = map(float,content[i].split()[:3])*info['dimension'] coords[:,i] = map(float,content[i].split()[:3])*info['size']
eulers[:,i] = map(float,content[i].split()[3:6]) eulers[:,i] = map(float,content[i].split()[3:6])
# -------------------------------------- switch according to task ---------------------------------- # -------------------------------------- switch according to task ----------------------------------
if options.config: if options.config: # write config file
file['output'].write('<microstructure>\n') file['output'].write('<microstructure>\n')
for i in xrange(info['grains']): for i in xrange(info['grains']):
file['output'].write('\n[Grain%s]\n'%(str(i+1).zfill(formatwidth)) + \ file['output'].write('\n[Grain%s]\n'%(str(i+1).zfill(formatwidth)) + \
@ -187,40 +174,56 @@ for file in files:
file['output'].write('\n[Grain%s]\n'%(str(i+1).zfill(formatwidth)) + \ file['output'].write('\n[Grain%s]\n'%(str(i+1).zfill(formatwidth)) + \
'(gauss)\tphi1 %g\tPhi %g\tphi2 %g\tscatter 0.0\tfraction 1.0\n'%(eulers[0,i],eulers[1,i],eulers[2,i])) '(gauss)\tphi1 %g\tPhi %g\tphi2 %g\tscatter 0.0\tfraction 1.0\n'%(eulers[0,i],eulers[1,i],eulers[2,i]))
else: else: # write geometry file
file['output'].write('%i\theader\n'%(len(new_header)) + ''.join(new_header)) N = info['grid'].prod()
shift = 0.5*info['size']/info['grid'] # shift by half of side length to center of element
N = info['resolution'].prod()
shift = 0.5*info['dimension']/info['resolution'] # shift by half of side length to center of element
undeformed = numpy.zeros((3,N),'d') undeformed = numpy.zeros((3,N),'d')
for i in xrange(N): for i in xrange(N):
undeformed[0,i] = info['dimension'][0]\ undeformed[0,i] = info['size'][0]\
* float(i % info['resolution'][0])\ * float(i % info['grid'][0])\
/float(info['resolution'][0]) /float(info['grid'][0])
undeformed[1,i] = info['dimension'][1]\ undeformed[1,i] = info['size'][1]\
* float(i//info['resolution'][0] % info['resolution'][1])\ * float(i//info['grid'][0] % info['grid'][1])\
/float(info['resolution'][1]) /float(info['grid'][1])
undeformed[2,i] = info['dimension'][2]\ undeformed[2,i] = info['size'][2]\
* float(i//info['resolution'][0]//info['resolution'][1] % info['resolution'][2])\ * float(i//info['grid'][0]//info['grid'][1] % info['grid'][2])\
/float(info['resolution'][2]) /float(info['grid'][2])
undeformed[:,i] += shift undeformed[:,i] += shift
indices = damask.core.math.periodicNearestNeighbor(\ indices = damask.core.math.periodicNearestNeighbor(\
info['dimension'],\ info['size'],\
numpy.eye(3),\ numpy.eye(3),\
undeformed,coords)//3**3 + 1 # floor division to kill periodic images undeformed,coords)//3**3 + 1 # floor division to kill periodic images
for n in xrange(info['resolution'][1:3].prod()): # loop over 2nd and 3rd dimension
file['output'].write({ True: ' ',
False:'\n'}[options.twoD].\
join(map(lambda x: str(x).rjust(formatwidth),\
indices[n*info['resolution'][0]:(n+1)*info['resolution'][0]]))+'\n')
missing = 0 missing = 0
for i in xrange(info['grains']): for i in xrange(info['grains']):
if i+1 not in indices: missing += 1 if i+1 not in indices: missing += 1
file['croak'].write({True:'all',False:'only'}[missing == 0] + ' %i grains mapped.\n'%(info['grains']-missing)) file['croak'].write({True:'all',False:'only'}[missing == 0] + ' %i grains mapped.\n'%(info['grains']-missing))
new_header.append("grid\ta %i\tb %i\tc %i\n"%(
info['grid'][0],
info['grid'][1],
info['grid'][2],))
new_header.append("size\tx %f\ty %f\tz %f\n"%(
info['size'][0],
info['size'][1],
info['size'][2],))
new_header.append("origin\tx %f\ty %f\tz %f\n"%(
info['origin'][0],
info['origin'][1],
info['origin'][2],))
new_header.append("microstructures\t%i\n"%(info['grains']-missing))
new_header.append("homogenization\t%i\n"%info['homogenization'])
file['output'].write('%i\theader\n'%(len(new_header)) + ''.join(new_header))
for n in xrange(info['grid'][1:3].prod()): # loop over 2nd and 3rd size
file['output'].write({ True: ' ',
False:'\n'}[options.twoD].\
join(map(lambda x: str(x).rjust(formatwidth),\
indices[n*info['grid'][0]:(n+1)*info['grid'][0]]))+'\n')
# ------------------------------------------ output finalization --------------------------------------- # ------------------------------------------ output finalization ---------------------------------------
if file['name'] != 'STDIN': if file['name'] != 'STDIN':

View File

@ -27,16 +27,16 @@ class extendedOption(Option):
# ----------------------- MAIN ------------------------------- # ----------------------- MAIN -------------------------------
identifiers = { identifiers = {
'resolution': ['a','b','c'], 'grid': ['a','b','c'],
'dimension': ['x','y','z'], 'size': ['x','y','z'],
'origin': ['x','y','z'], 'origin': ['x','y','z'],
} }
mappings = { mappings = {
'resolution': lambda x: int(x), 'grid': lambda x: int(x),
'dimension': lambda x: float(x), 'size': lambda x: float(x),
'origin': lambda x: float(x), 'origin': lambda x: float(x),
'homogenization': lambda x: int(x), 'homogenization': lambda x: int(x),
'maxmicrostructure': lambda x: int(x), 'microstructures': lambda x: int(x),
} }
@ -51,19 +51,26 @@ compress geometry files with ranges "a to b" and/or multiples "n of x".
files = [] files = []
if filenames == []: if filenames == []:
files.append({'name':'STDIN', 'input':sys.stdin, 'output':sys.stdout}) files.append({'name':'STDIN',
'input':sys.stdin,
'output':sys.stdout,
'croak':sys.stderr,
})
else: else:
for name in filenames: for name in filenames:
if os.path.exists(name): if os.path.exists(name):
files.append({'name':name, 'input':open(name), 'output':open(name+'_tmp','w')}) files.append({'name':name,
'input':open(name),
'output':open(name+'_tmp','w'),
'croak':sys.stdout,
})
# ------------------------------------------ loop over input files --------------------------------------- # ------------------------------------------ loop over input files ---------------------------------------
for file in files: for file in files:
if file['name'] != 'STDIN': print file['name'] if file['name'] != 'STDIN': file['croak'].write(file['name']+'\n')
# get labels by either read the first row, or - if keyword header is present - the last line of the header # get labels by either read the first row, or - if keyword header is present - the last line of the header
firstline = file['input'].readline() firstline = file['input'].readline()
m = re.search('(\d+)\s*head', firstline.lower()) m = re.search('(\d+)\s*head', firstline.lower())
if m: if m:
@ -76,16 +83,18 @@ for file in files:
content = file['input'].readlines() content = file['input'].readlines()
file['input'].close() file['input'].close()
info = {'resolution': [0,0,0], info = {'grid': [0,0,0],
'dimension': [0.0,0.0,0.0], 'size': [0.0,0.0,0.0],
'origin': [0.0,0.0,0.0], 'origin': [0.0,0.0,0.0],
'homogenization': 1, 'homogenization': 1,
'maxmicrostructure': 0, 'microstructures': 0,
} }
new_header = [] new_header = []
for header in headers: for header in headers:
headitems = map(str.lower,header.split()) headitems = map(str.lower,header.split())
if headitems[0] == 'resolution': headitems[0] = 'grid'
if headitems[0] == 'dimension': headitems[0] = 'size'
if headitems[0] in mappings.keys(): if headitems[0] in mappings.keys():
if headitems[0] in identifiers.keys(): if headitems[0] in identifiers.keys():
for i in xrange(len(identifiers[headitems[0]])): for i in xrange(len(identifiers[headitems[0]])):
@ -94,33 +103,34 @@ for file in files:
else: else:
info[headitems[0]] = mappings[headitems[0]](headitems[1]) info[headitems[0]] = mappings[headitems[0]](headitems[1])
if info['resolution'] == [0,0,0]: if info['grid'] == [0,0,0]:
sys.stderr.write('no resolution info found.\n') file['croak'].write('no grid info found.\n')
continue continue
if info['dimension'] == [0.0,0.0,0.0]: if info['size'] == [0.0,0.0,0.0]:
sys.stderr.write('no dimension info found.\n') file['croak'].write('no size info found.\n')
continue continue
if file['name'] != 'STDIN': file['croak'].write('grid a b c: %s\n'%(' x '.join(map(str,info['grid']))) + \
print 'resolution: %s'%(' x '.join(map(str,info['resolution']))) 'size x y z: %s\n'%(' x '.join(map(str,info['size']))) + \
print 'dimension: %s'%(' x '.join(map(str,info['dimension']))) 'origin x y z: %s\n'%(' : '.join(map(str,info['origin']))) + \
print 'origin: %s'%(' : '.join(map(str,info['origin']))) 'homogenization: %i\n'%info['homogenization'] + \
'microstructures: %i\n'%info['microstructures'])
new_header.append("resolution\ta %i\tb %i\tc %i\n"%( new_header.append("grid\ta %i\tb %i\tc %i\n"%(
info['resolution'][0], info['grid'][0],
info['resolution'][1], info['grid'][1],
info['resolution'][2],)) info['grid'][2],))
new_header.append("dimension\tx %f\ty %f\tz %f\n"%( new_header.append("size\tx %f\ty %f\tz %f\n"%(
info['dimension'][0], info['size'][0],
info['dimension'][1], info['size'][1],
info['dimension'][2])) info['size'][2]))
new_header.append("origin\tx %f\ty %f\tz %f\n"%( new_header.append("origin\tx %f\ty %f\tz %f\n"%(
info['origin'][0], info['origin'][0],
info['origin'][1], info['origin'][1],
info['origin'][2])) info['origin'][2]))
new_header.append("homogenization\t%i\n"%info['homogenization']) new_header.append("homogenization\t%i\n"%info['homogenization'])
if info['maxmicrostructure'] > 0: if info['microstructures'] > 0:
new_header.append("maxMicrostructure\t%i\n"%info['maxmicrostructure']) new_header.append("microstructures\t%i\n"%info['microstructures'])
# ------------------------------------------ assemble header --------------------------------------- # ------------------------------------------ assemble header ---------------------------------------

View File

@ -27,16 +27,17 @@ class extendedOption(Option):
# ----------------------- MAIN ------------------------------- # ----------------------- MAIN -------------------------------
identifiers = { identifiers = {
'resolution': ['a','b','c'], 'grid': ['a','b','c'],
'dimension': ['x','y','z'], 'size': ['x','y','z'],
'origin': ['x','y','z'], 'origin': ['x','y','z'],
} }
mappings = { mappings = {
'resolution': lambda x: int(x), 'grid': lambda x: int(x),
'dimension': lambda x: float(x), 'size': lambda x: float(x),
'origin': lambda x: float(x), 'origin': lambda x: float(x),
'homogenization': lambda x: int(x), 'homogenization': lambda x: int(x),
'maxmicrostructure': lambda x: int(x), 'microstructures': lambda x: int(x),
} }
@ -46,26 +47,33 @@ Unpack geometry files containing ranges "a to b" and/or "n of x" multiples (excl
) )
parser.add_option('-2', '--twodimensional', dest='twoD', action='store_true', \ parser.add_option('-2', '--twodimensional', dest='twoD', action='store_true', \
help='output geom file with two-dimensional data arrangement') help='output geom file with two-dimensional data arrangement [%default]')
parser.set_defaults(twoD = False) parser.set_defaults(twoD = False)
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
# ------------------------------------------ setup file handles --------------------------------------- # ------------------------------------------ setup file handles ---------------------------------------
files = [] files = []
if filenames == []: if filenames == []:
files.append({'name':'STDIN', 'input':sys.stdin, 'output':sys.stdout}) files.append({'name':'STDIN',
'input':sys.stdin,
'output':sys.stdout,
'croak':sys.stderr,
})
else: else:
for name in filenames: for name in filenames:
if os.path.exists(name): if os.path.exists(name):
files.append({'name':name, 'input':open(name), 'output':open(name+'_tmp','w')}) files.append({'name':name,
'input':open(name),
'output':open(name+'_tmp','w'),
'croak':sys.stdout,
})
# ------------------------------------------ loop over input files --------------------------------------- # ------------------------------------------ loop over input files ---------------------------------------
for file in files: for file in files:
if file['name'] != 'STDIN': print file['name'] if file['name'] != 'STDIN': file['croak'].write(file['name']+'\n')
# get labels by either read the first row, or - if keyword header is present - the last line of the header # get labels by either read the first row, or - if keyword header is present - the last line of the header
@ -81,16 +89,18 @@ for file in files:
content = file['input'].readlines() content = file['input'].readlines()
file['input'].close() file['input'].close()
info = {'resolution': [0,0,0], info = {'grid': [0,0,0],
'dimension': [0.0,0.0,0.0], 'size': [0.0,0.0,0.0],
'origin': [0.0,0.0,0.0], 'origin': [0.0,0.0,0.0],
'homogenization': 1, 'homogenization': 1,
'maxmicrostructure': 0, 'microstructures': 0,
} }
new_header = [] new_header = []
for header in headers: for header in headers:
headitems = map(str.lower,header.split()) headitems = map(str.lower,header.split())
if headitems[0] == 'resolution': headitems[0] = 'grid'
if headitems[0] == 'dimension': headitems[0] = 'size'
if headitems[0] in mappings.keys(): if headitems[0] in mappings.keys():
if headitems[0] in identifiers.keys(): if headitems[0] in identifiers.keys():
for i in xrange(len(identifiers[headitems[0]])): for i in xrange(len(identifiers[headitems[0]])):
@ -99,39 +109,40 @@ for file in files:
else: else:
info[headitems[0]] = mappings[headitems[0]](headitems[1]) info[headitems[0]] = mappings[headitems[0]](headitems[1])
if info['resolution'] == [0,0,0]: if info['grid'] == [0,0,0]:
sys.stderr.write('no resolution info found.\n') file['croak'].write('no grid info found.\n')
continue continue
if info['dimension'] == [0.0,0.0,0.0]: if info['size'] == [0.0,0.0,0.0]:
sys.stderr.write('no dimension info found.\n') file['croak'].write('no size info found.\n')
continue continue
format = {True: info['resolution'][0], format = {True: info['grid'][0],
False: 1}[options.twoD] False: 1}[options.twoD]
if file['name'] != 'STDIN': file['croak'].write('grid a b c: %s\n'%(' x '.join(map(str,info['grid']))) + \
print 'resolution: %s'%(' x '.join(map(str,info['resolution']))) 'size x y z: %s\n'%(' x '.join(map(str,info['size']))) + \
print 'dimension: %s'%(' x '.join(map(str,info['dimension']))) 'origin x y z: %s\n'%(' : '.join(map(str,info['origin']))) + \
print 'origin: %s'%(' : '.join(map(str,info['origin']))) 'homogenization: %i\n'%info['homogenization'] + \
'microstructures: %i\n'%info['microstructures'])
new_header.append("resolution\ta %i\tb %i\tc %i\n"%( new_header.append("grid\ta %i\tb %i\tc %i\n"%(
info['resolution'][0], info['grid'][0],
info['resolution'][1], info['grid'][1],
info['resolution'][2],)) info['grid'][2],))
new_header.append("dimension\tx %f\ty %f\tz %f\n"%( new_header.append("size\tx %f\ty %f\tz %f\n"%(
info['dimension'][0], info['size'][0],
info['dimension'][1], info['size'][1],
info['dimension'][2])) info['size'][2]))
new_header.append("origin\tx %f\ty %f\tz %f\n"%( new_header.append("origin\tx %f\ty %f\tz %f\n"%(
info['origin'][0], info['origin'][0],
info['origin'][1], info['origin'][1],
info['origin'][2])) info['origin'][2]))
new_header.append("homogenization\t%i\n"%info['homogenization']) new_header.append("homogenization\t%i\n"%info['homogenization'])
if info['maxmicrostructure'] > 0: if info['microstructures'] > 0:
new_header.append("maxMicrostructure\t%i\n"%info['maxmicrostructure']) new_header.append("microstructures\t%i\n"%info['microstructures'])
digits = 1+int(math.log10(int(info['maxmicrostructure']))) digits = 1+int(math.log10(int(info['microstructures'])))
else: else:
digits = 1+int(math.log10(int(info['resolution'][0]*info['resolution'][1]*info['resolution'][2]))) digits = 1+int(math.log10(int(info['grid'][0]*info['grid'][1]*info['grid'][2])))
# ------------------------------------------ assemble header --------------------------------------- # ------------------------------------------ assemble header ---------------------------------------

View File

@ -27,12 +27,10 @@ class extendedOption(Option):
# ----------------------- MAIN ------------------------------- # ----------------------- MAIN -------------------------------
identifiers = { identifiers = {
'resolution': ['a','b','c'], 'grid': ['a','b','c'],
'dimension': ['x','y','z'],
} }
mappings = { mappings = {
'resolution': lambda x: int(x), 'grid': lambda x: int(x),
'dimension': lambda x: float(x),
} }
@ -44,18 +42,18 @@ Reports positions with random crystal orientations in seeds file format to STDOU
parser.add_option('-N', dest='N', type='int', \ parser.add_option('-N', dest='N', type='int', \
help='number of seed points to distribute [%default]') help='number of seed points to distribute [%default]')
parser.add_option('-r','--resolution', dest='res', type='int', nargs=3, \ parser.add_option('-g','--grid', dest='grid', type='int', nargs=3, \
help='Min Fourier points in x, y, z %default') help='min a,b,c grid of hexahedral box %default')
parser.add_option('-s', '--rnd', dest='randomSeed', type='int', \ parser.add_option('-r', '--rnd', dest='randomSeed', type='int', \
help='seed of random number generator [%default]') help='seed of random number generator [%default]')
parser.set_defaults(randomSeed = 0) parser.set_defaults(randomSeed = 0)
parser.set_defaults(res = [16,16,16]) parser.set_defaults(grid = [16,16,16])
parser.set_defaults(N = 20) parser.set_defaults(N = 20)
(options, extras) = parser.parse_args() (options, extras) = parser.parse_args()
Npoints = options.res[0]*options.res[1]*options.res[2] Npoints = options.grid[0]*options.grid[1]*options.grid[2]
if options.N > Npoints: if options.N > Npoints:
sys.stderr.write('Warning: more seeds than grid points at minimum resolution.\n') sys.stderr.write('Warning: more seeds than grid points at minimum resolution.\n')
options.N = Npoints options.N = Npoints
@ -69,14 +67,17 @@ grainEuler[1,:] = numpy.arccos(2*grainEuler[1,:]-1)*180.0/math.pi
grainEuler[2,:] *= 360.0 grainEuler[2,:] *= 360.0
seedpoint = numpy.random.permutation(Npoints)[:options.N] seedpoint = numpy.random.permutation(Npoints)[:options.N]
seeds[0,:] = (numpy.mod(seedpoint ,options.res[0])+numpy.random.random())/options.res[0] seeds[0,:]=(numpy.mod(seedpoint ,options.grid[0])\
seeds[1,:] = (numpy.mod(seedpoint// options.res[0] ,options.res[1])+numpy.random.random())/options.res[1] +numpy.random.random())/options.grid[0]
seeds[2,:] = (numpy.mod(seedpoint//(options.res[1]*options.res[0]),options.res[2])+numpy.random.random())/options.res[2] seeds[1,:]=(numpy.mod(seedpoint// options.grid[0],options.grid[1])\
+numpy.random.random())/options.grid[1]
seeds[2,:]=(numpy.mod(seedpoint//(options.grid[1]*options.grid[0]),options.grid[2])\
+numpy.random.random())/options.grid[2]
print "4\theader" print "4\theader"
print "resolution\ta %i\tb %i\tc %i"%(options.res[0],options.res[1],options.res[2],) print "grid\ta %i\tb %i\tc %i"%(options.grid[0],options.grid[1],options.grid[2],)
print "grains\t%i"%options.N print "grains\t%i"%options.N
print "randomSeed\t%f"%(options.randomSeed) print "randomSeed\t%i"%(options.randomSeed)
print "x\ty\tz\tphi1\tPhi\tphi2" print "x\ty\tz\tphi1\tPhi\tphi2"
numpy.savetxt(sys.stdout,numpy.transpose(numpy.concatenate((seeds,grainEuler),axis = 0)),fmt='%10.6f',delimiter='\t') numpy.savetxt(sys.stdout,numpy.transpose(numpy.concatenate((seeds,grainEuler),axis = 0)),fmt='%10.6f',delimiter='\t')