Python FortranFile.read_record方法代码示例

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
def read():
	f=FortranFile(BinaryData,'r')
	n_vals=f.read_record('i4,i4,i4,i4,i4,(12,)i4')
	d_vals=f.read_reals('f4')
	xyz_vals=f.read_record('(500,)f4,(500,)f4,(500,)f4')
	a_vals=f.read_reals('i8')
	ndr_vals=f.read_ints('i4')
	ag_vals=f.read_record('(3,)i8')
	f.close()
	return n_vals,d_vals,xyz_vals,a_vals,ndr_vals,ag_vals

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
def reorder_uXu(ext,formatted=False):
  try:
    print "----------\n  {0}  \n----------".format(ext)
    
    if formatted:
	f_uXu_in = open(seedname+"."+ext, 'r')
	f_uXu_out = open(seednameGW+"."+ext, 'w')
	header=f_uXu_in.readline()
	f_uXu_out.write(header)
	nbnd,NK,nnb=np.array(f_uXu_in.readline().split(),dtype=int)
	f_uXu_out.write("  ".join(str(x) for x in [NBND,NK,nnb])+"\n")
    else:
	f_uXu_in = FortranFile(seedname+"."+ext, 'r')
	f_uXu_out = FortranFile(seednameGW+"."+ext, 'w')
	header=f_uXu_in.read_record(dtype='c')
        f_uXu_out.write_record(header)
	nbnd,NK,nnb=np.array(f_uXu_in.read_record(dtype=np.int32))
	f_uXu_out.write_record(np.array([NBND,NK,nnb],dtype=np.int32))
    
    assert nbnd==nbndDFT
    print nbnd,NK,nnb
    
    if formatted:
	uXu=np.loadtxt(f_uXu_in).reshape(-1)
	start=0
	length=nbnd*nbnd

    for ik in xrange(NKPT):
	for ib2 in xrange(nnb):
	    for ib1 in xrange(nnb):
		if formatted:
		    A=uXu[start:start+length]
		    start+=length
		else:
		    A=f_uXu_in.read_record(dtype=np.complex)
		A=(A.reshape(nbnd,nbnd,order='F')[BANDSORT[KPNB[ik][ib2]],:][:,BANDSORT[KPNB[ik][ib1]]]+
			np.einsum('ln,lm,l->nm',MMN[ik][ib2].conj(),MMN[ik][ib1],eigenDE[ik]) ).reshape(-1,order='F')
		if formatted:
		    f_uXu_out.write("".join("{0:26.16e}  {1:26.16f}\n".format(x.real,x.imag) for x in A))
		else:
	    	    f_uXu_out.write_record(A)
    f_uXu_out.close()
    f_uXu_in.close()
    print "----------\n {0} OK  \n----------\n".format(ext)
  except IOError as err:
    print "WARNING: {0}.{1} not written : ".format(seednameGW,ext),err 

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
def test_fortranfiles_mixed_record():
    filename = path.join(DATA_PATH, "fortran-mixed.dat")
    f = FortranFile(filename, 'r', '<u4')
    record = f.read_record('<i4,<f4,<i8,(2)<f8')
    assert_equal(record['f0'][0], 1)
    assert_allclose(record['f1'][0], 2.3)
    assert_equal(record['f2'][0], 4)
    assert_allclose(record['f3'][0], [5.6, 7.8])

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
	def read_stalker_data(self,
		datadir='/data',
		procdir='proc0'):
			
		f = FortranFile(datadir+'/'+procdir+'/particles_stalker.dat','r')
		file_not_ended=True
		firstline = f.read_record([('a','f8'),('b','i4')])
		firstline = f.read_record([('a','f8'),('b','i4')])
		snaptime = []
		nrshot = []
		partdata = []
		iddata = []
		while (True):
			try:
				temp = f.read_record([('snaptime','f8'),('nrshot','i4')])
				snaptime.append(temp['snaptime'])
				nrshot.append(temp['nrshot'])
			except TypeError:
				break
			except ValueError:
				break
			if (nrshot[-1] > 0):
				temp = f.read_record([('iddata','i4')])
				iddata.append(temp['iddata'])
				temp = f.read_record([('partdata','f8')])
				temp = np.reshape(temp,(-1,self.ndims_stalk))
				partdata.append(temp['partdata'])
			else:
				iddata.append([])
				partdata.append([])
		
		partdata=np.asarray(partdata)
		iddata = np.asarray(iddata)
		snaptime = np.asarray(snaptime)
		nrshot = np.array(nrshot)
		oldshape = np.shape(partdata)
		
		partdata = np.reshape(partdata,oldshape)
		
		setattr(self,procdir+'_snaptime',snaptime)
		setattr(self,procdir+'_iddata',iddata)
		setattr(self,procdir+'_partdata',partdata)
		setattr(self,procdir+'_nrshot',nrshot)

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
def test_fortranfiles_read():
    for filename in iglob(path.join(DATA_PATH, "fortran-*-*x*x*.dat")):
        m = re.search(r'fortran-([^-]+)-(\d+)x(\d+)x(\d+).dat', filename, re.I)
        if not m:
            raise RuntimeError("Couldn't match %s filename to regex" % filename)

        dims = (int(m.group(2)), int(m.group(3)), int(m.group(4)))

        dtype = m.group(1).replace('s', '<')

        f = FortranFile(filename, 'r', '<u4')
        data = f.read_record(dtype=dtype).reshape(dims, order='F')
        f.close()

        expected = np.arange(np.prod(dims)).reshape(dims).astype(dtype)
        assert_equal(data, expected)

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
def test_fortranfiles_read():
    for filename in iglob(path.join(DATA_PATH, "fortran-*-*x*x*.dat")):
        m = re.search('fortran-([^-]+)-(\d+)x(\d+)x(\d+).dat', filename, re.I)
        if not m:
            raise RuntimeError("Couldn't match %s filename to regex" % filename)
        dims = (int(m.group(2)), int(m.group(3)), int(m.group(4)))

        f = FortranFile(filename, 'r', '<u4')
        data = f.read_record(dtype=m.group(1)).reshape(dims)
        f.close()

        counter = 0
        for k in range(dims[2]):
            for j in range(dims[1]):
                for i in range(dims[0]):
                    assert_equal(counter, data[i,j,k])
                    counter += 1

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
    def read(self, datadir='data', proc=-1, quiet=False,
             trim=False):
        """
        Read the grid data from the pencil code simulation.
        If proc < 0, then load all data and assemble.
        Otherwise, load grid from specified processor.

        call signature:

        read(self, file_name='time_series.dat', datadir='data',
             double=0, quiet=0, comment_char='#')

        Keyword arguments:

        *datadir*:
          Directory where the data is stored.

        *proc*
          Processor to be read. If proc is -1, then read the 'global'
          grid. If proc is >=0, then read the grid.dat in the
          corresponding processor directory.

        *quiet*
          Flag for switching of output.

        *trim*
          Cuts off the ghost points.
        """

        import numpy as np
        import os
        from scipy.io import FortranFile
        import pencilnew.read as read

        datadir = os.path.expanduser(datadir)
        dim = read.dim(datadir, proc)
        if dim.precision == 'D':
            precision = 'd'
        else:
            precision = 'f'

        if proc < 0:
            proc_dirs = list(filter(lambda string: string.startswith('proc'), os.listdir(datadir)))
        else:
            proc_dirs = ['proc' + str(proc)]

        # Define the global arrays.
        x = np.zeros(dim.mx, dtype=precision)
        y = np.zeros(dim.my, dtype=precision)
        z = np.zeros(dim.mz, dtype=precision)
        dx_1 = np.zeros(dim.mx, dtype=precision)
        dy_1 = np.zeros(dim.my, dtype=precision)
        dz_1 = np.zeros(dim.mz, dtype=precision)
        dx_tilde = np.zeros(dim.mx, dtype=precision)
        dy_tilde = np.zeros(dim.my, dtype=precision)
        dz_tilde = np.zeros(dim.mz, dtype=precision)

        for directory in proc_dirs:
            proc = int(directory[4:])
            procdim = read.dim(datadir, proc)
            if not quiet:
                print("reading grid data from processor {0} of {1} ...".format(proc, len(proc_dirs)))

            mxloc = procdim.mx
            myloc = procdim.my
            mzloc = procdim.mz

            # Read the grid data.
            file_name = os.path.join(datadir, directory, 'grid.dat')
            infile = FortranFile(file_name, 'r')
            grid_raw = infile.read_record(dtype=precision)
            dx, dy, dz = tuple(infile.read_record(dtype=precision))
            Lx, Ly, Lz = tuple(infile.read_record(dtype=precision))
            dx_1_raw = infile.read_record(dtype=precision)
            dx_tilde_raw = infile.read_record(dtype=precision)
            infile.close()

            # Reshape the arrays.
            t = grid_raw[0]
            x_loc = grid_raw[1:mxloc+1]
            y_loc = grid_raw[mxloc+1:mxloc+myloc+1]
            z_loc = grid_raw[mxloc+myloc+1:mxloc+myloc+mzloc+1]
            dx_1_loc = dx_1_raw[0:mxloc]
            dy_1_loc = dx_1_raw[mxloc:mxloc+myloc]
            dz_1_loc = dx_1_raw[mxloc+myloc:mxloc+myloc+mzloc]
            dx_tilde_loc = dx_tilde_raw[0:mxloc]
            dy_tilde_loc = dx_tilde_raw[mxloc:mxloc+myloc]
            dz_tilde_loc = dx_tilde_raw[mxloc+myloc:mxloc+myloc+mzloc]

            if len(proc_dirs) > 1:
                if procdim.ipx == 0:
                    i0x = 0
                    i1x = i0x + procdim.mx
                    i0x_loc = 0
                    i1x_loc = procdim.mx
                else:
                    i0x = procdim.ipx*procdim.nx + procdim.nghostx
                    i1x = i0x + procdim.mx - procdim.nghostx
                    i0x_loc = procdim.nghostx
                    i1x_loc = procdim.mx
#.........这里部分代码省略.........

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
def main(sourcedir, sourcefile, savedir, saveheader, moving=0, incr=1, imin=0, imax=-1 ):
    """
    imin --> minimum iteration above which files will be written
    incr --> iteration interval to write files at
    """

    # ogdir = os.getcwd()
    # os.chdir(savedir)

    MakeOutputDir(savedir)

    #link in source file
    if not os.path.isfile('{}/{}'.format(savedir, sourcefile)):
        # cmd('ln -s {}/{} {}'.format(sourcedir, sourcefile, sourcefile))
        cmd('ln -s {}/{} {}/{}'.format(sourcedir, sourcefile, savedir, sourcefile))


    #GET GRID/FILE DIMENSIONS
    f = FortranFile('{}/{}'.format(savedir, sourcefile), 'r')
    dims = f.read_ints('uint32')[2:6]
    #Number of iterations, Number of points in each direction, number of parameters, needed to read whole file
    nj, nk, nl, nq = dims

    #FORMATS IN BC201 FILE
    bc201formats = [
                        ('ints', '<i4', (1,7)), #IGRID,ISTEP,NJ,NK,NL,NQ,NQC
                        ('tvr', 'float64'),     #TVREF
                        ('dtr', 'float64'),     #DTVREF
                        ('xyz', 'float64', (3,nj,nk,nl)), #XYZ points
                        ('q', 'float64', (nq,nj,nk,nl)),  #Q variables (nq variables, njXnkXnl values for each)
                        ('iblank', 'int32', (nj,nk,nl)),  #IBLANK status of each grid point
                    ]

    #########################
    #Save pressure coeff history at certian locations
    #! I want to write out a subset through time (every iteration!)
    ports = [8,13,30,45]
    out = open(savedir+"/cp_history.dat","w")
    out.write("ITER")
    for p in ports:
        out.write(" {}_x {}_cp".format(p,p))
    out.write("\n")
    #############################


    #RESTART FILE AND READ EACH RECORD
    f = FortranFile('{}/{}'.format(savedir, sourcefile), 'r')
    keep_reading = True
    irecord = -1
    while (keep_reading is True):

        irecord += 1

        #READ CURRENT RECORD
        b = f.read_record(bc201formats)



        #GET CURRENT ITERATION
        istep = b['ints'][0][0][1]

        # print(istep)


        #DONT WRITE UNDESIRED FILES (SKIP)
        #only write files at specified interval
        if (istep % incr != 0):
            continue
        #Don't write files below start iter
        if (istep < imin):
            continue
        #Stop reading when max desired iteration is reached
        if istep > imax - incr:
            keep_reading = False




        #WRITE GRID POINTS TO PLOT3D FORMATTED FILE

        if moving:
            #Grid is moving, write to file at each iteration
            savename = '{}/{}.x.{}'.format(savedir, saveheader, istep)
            WriteGrid(savename, b, nj, nk, nl)
        elif irecord == 0:
           #Grid is stationary, write to file only once
            savename = '{}/{}.x'.format(savedir, saveheader)
            WriteGrid(savename, b, nj, nk, nl)


        #CALCULATE CP FOR EACH POINT
        cps = np.zeros((nj,nk,nl))
        # print(cps.shape)
        for j in range(nj):
            for k in range(nk):
                for l in range(nl):
                    #print(j,k,l,b['q'].shape)
                    q = np.zeros(nq)
                    for iq in range(nq):
                        q[iq] = b['q'][0][iq][j][k][l]
#.........这里部分代码省略.........

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]
if calcSPN:
  try:
    print "----------\n SPN  \n---------\n"

    if SPNformatted:
	f_spn_in =  open(seedname+".spn", 'r')
	f_spn_out = open(seednameGW+".spn", 'w')
	header=f_spn_in.readline()
        f_spn_out.write(header)
        nbnd,NK=np.array(f_spn_in.readline().split(),dtype=np.int32)
	f_spn_out.write("  ".join(str(x) for x in (NBND,NKPT) ) )
    else:
	f_spn_in = FortranFile(seedname+".spn", 'r')
	f_spn_out = FortranFile(seednameGW+".spn", 'w')
	header=f_spn_in.read_record(dtype='c')
        f_spn_out.write_record(header)
        nbnd,NK=f_spn_in.read_record(dtype=np.int32)
	f_spn_out.write_record(np.array([NBND,NKPT],dtype=np.int32))

    print "".join(header)
    assert nbnd==nbndDFT


    indm,indn=np.tril_indices(nbnd)
    indmQP,indnQP=np.tril_indices(NBND)

    if SPNformatted:
	SPN=np.loadtxt(f_spn_in).reshape(-1)
	start=0
	length=(3*nbnd*(nbnd+1))/2

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]

#.........这里部分代码省略.........

        # Set up the global array.
        if not run2D:
            f = np.zeros((total_vars, dim.mz, dim.my, dim.mx),
                         dtype=precision)
        else:
            if dim.ny == 1:
                f = np.zeros((total_vars, dim.mz, dim.mx), dtype=precision)
            else:
                f = np.zeros((total_vars, dim.my, dim.mx), dtype=precision)

        x = np.zeros(dim.mx, dtype=precision)
        y = np.zeros(dim.my, dtype=precision)
        z = np.zeros(dim.mz, dtype=precision)

        for directory in proc_dirs:
            proc = int(directory[4:])
            if var_file[0:2].lower() == 'og':
                procdim = read.ogdim(datadir, proc)
            else:
                procdim = read.dim(datadir, proc)
            if not quiet:
                print("Reading data from processor {0} of {1} ...".format( \
                      proc, len(proc_dirs)))

            mxloc = procdim.mx
            myloc = procdim.my
            mzloc = procdim.mz

            # Read the data.
            file_name = os.path.join(datadir, directory, var_file)
            infile = FortranFile(file_name)
            if not run2D:
                f_loc = infile.read_record(dtype=precision)
                f_loc = f_loc.reshape((-1, mzloc, myloc, mxloc))
            else:
                if dim.ny == 1:
                    f_loc = infile.read_record(dtype=precision)
                    f_loc = f_loc.reshape((-1, mzloc, mxloc))
                else:
                    f_loc = infile.read_record(dtype=precision)
                    f_loc = f_loc.reshape((-1, myloc, mxloc))
            raw_etc = infile.read_record(precision)
            infile.close()

            t = raw_etc[0]
            x_loc = raw_etc[1:mxloc+1]
            y_loc = raw_etc[mxloc+1:mxloc+myloc+1]
            z_loc = raw_etc[mxloc+myloc+1:mxloc+myloc+mzloc+1]
            if param.lshear:
                shear_offset = 1
                deltay = raw_etc[-1]
            else:
                shear_offset = 0

            dx = raw_etc[-3-shear_offset]
            dy = raw_etc[-2-shear_offset]
            dz = raw_etc[-1-shear_offset]

            if len(proc_dirs) > 1:
                # Calculate where the local processor will go in
                # the global array.
                #
                # Don't overwrite ghost zones of processor to the left (and
                # accordingly in y and z direction -- makes a difference on the
                # diagonals)

# 需要导入模块: from scipy.io import FortranFile [as 别名]
# 或者: from scipy.io.FortranFile import read_record [as 别名]

#.........这里部分代码省略.........
            try:
                field_list.remove('position')
            except:
                pass

        # Initialize the extensions list.
        if extension:
            if isinstance(extension, list):
                extension_list = extension
            else:
                extension_list = [extension]
        else:
            # Find the existing extensions.
            extension_list = []
            for file_name in os.listdir(slice_dir):
                if file_name[:6] == 'slice_':
                    extension_list.append(file_name.split('.')[1])
            # Remove duplicates.
            extension_list = list(set(extension_list))
            try:
                extension_list.remove('dat')
            except:
                pass

        class Foo(object):
            pass

        for extension in extension_list:
            if verbose:
                print('Extension: ' + str(extension))
            # This one will store the data.
            ext_object = Foo()

            for field in field_list:
                if verbose:
                    print('  -> Field: ' + str(field))
                # Compose the file name according to field and extension.
                datadir = os.path.expanduser(datadir)
                if proc < 0:
                    file_name = os.path.join(datadir, 'slice_'+field+'.'+extension)
                else:
                    file_name = os.path.join(datadir, 'proc{0}'.format(proc),
                                             'slice_'+field+'.'+extension)

                dim = read.dim(datadir, proc)
                if dim.precision == 'D':
                    precision = 'd'
                else:
                    precision = 'f'

                # Set up slice plane.
                if extension == 'xy' or extension == 'Xy' or  extension == 'xy2':
                    hsize = dim.nx
                    vsize = dim.ny
                if extension == 'xz':
                    hsize = dim.nx
                    vsize = dim.nz
                if extension == 'yz':
                    hsize = dim.ny
                    vsize = dim.nz

                try:
                    infile = FortranFile(file_name)
                except:
                    continue

                islice = 0
                self.t = np.zeros(1, dtype=precision)
                self.t = [0]
                slice_series = [0]

                while True:
                    try:
                        if verbose:
                            print('  -> Reading... ')
                        raw_data = infile.read_record(dtype=precision)
                    except ValueError:
                        break
                    except TypeError:
                        break

                    if old_file:
                        self.t.append(list(raw_data[-1]))
                        slice_series.extend(list(raw_data[:-1]))
                    else:
                        self.t.append(list(raw_data[-2:-1]))
                        slice_series.extend(list(raw_data[:-2]))
                    islice += 1
                    if verbose:
                        print('  -> Done')

                # Reshape and remove first entry.
                if verbose:
                    print('Reshaping array')
                self.t = np.array(self.t[1:], dtype=precision)[:, 0]
                slice_series = np.array(slice_series, dtype=precision)
                slice_series = slice_series[1:].reshape(islice, vsize, hsize)
                setattr(ext_object, field, slice_series)

            setattr(self, extension, ext_object)