source: inundation/ga/storm_surge/pyvolution/data_manager.py @ 668

Last change on this file since 668 was 620, checked in by ole, 20 years ago

Added functionality for storing momentum (x and y) as well as water levels in sww files.
File size: 24.3 KB
Line 
1"""Functions to store and retrieve data for the Shallow Water Wave equation.
2There are two kinds of data
3
4  1: Constant data: Vertex coordinates and field values. Stored once
5  2: Variable data: Conserved quantities. Stored once per timestep.
6
7All data is assumed to reside at vertex locations. 
8
9"""
10
11from Numeric import concatenate       
12
13
14def make_filename(s):
15    """Transform argument string into a suitable filename
16    """
17
18    s = s.strip()
19    s = s.replace(' ', '_')
20    s = s.replace('(', '')
21    s = s.replace(')', '')
22    s = s.replace('__', '_')               
23       
24    return s
25
26
27def check_dir(path, verbose=None):
28    """Check that specified path exists.
29    If path does not exist it will be created if possible   
30
31    USAGE:
32       checkdir(path, verbose):
33
34    ARGUMENTS:
35        path -- Directory
36        verbose -- Flag verbose output (default: None)
37
38    RETURN VALUE:
39        Verified path including trailing separator
40       
41    """
42
43    import os, sys
44    import os.path
45
46    if sys.platform in ['nt', 'dos', 'win32', 'what else?']:
47        unix = 0
48    else:
49        unix = 1
50
51
52    if path[-1] != os.sep: 
53        path = path + os.sep  # Add separator for directories
54
55    path = os.path.expanduser(path) # Expand ~ or ~user in pathname
56    if not (os.access(path,os.R_OK and os.W_OK) or path == ''):
57        try:
58            exitcode=os.mkdir(path)
59
60            # Change access rights if possible
61            #
62            if unix:
63                exitcode=os.system('chmod 775 '+path)
64            else:
65                pass  # FIXME: What about acces rights under Windows?
66           
67            if verbose: print 'MESSAGE: Directory', path, 'created.'
68         
69        except:
70            print 'WARNING: Directory', path, 'could not be created.'
71            if unix:
72                path = '/tmp/'
73            else:
74                path = 'C:'
75               
76            print 'Using directory %s instead' %path
77
78    return(path)
79
80   
81
82def del_dir(path):
83    """Recursively delete directory path and all its contents
84    """
85
86    import os
87
88    if os.path.isdir(path):
89        for file in os.listdir(path):
90            X = os.path.join(path, file)
91
92
93            if os.path.isdir(X) and not os.path.islink(X):
94                del_dir(X) 
95            else:
96                try:
97                    os.remove(X)
98                except:
99                    print "Could not remove file %s" %X
100
101        os.rmdir(path)
102
103
104       
105def create_filename(filename, format, size, time=None):
106
107    import os
108    from config import data_dir
109
110           
111    FN = check_dir(data_dir) + filename + '_size%d' %size
112
113    if time is not None:
114        FN += '_time%.2f' %time
115       
116    FN += '.' + format
117    return FN
118   
119
120def get_files(filename, format, size):
121    """Get all file (names) with gven name, size and format
122    """
123
124    import glob
125
126    import os
127    from config import data_dir
128
129           
130    dir = check_dir(data_dir)
131
132    pattern = dir + os.sep + filename + '_size=%d*.%s' %(size, format)
133    return glob.glob(pattern)
134
135
136
137#Generic class for storing output to e.g. visualisation or checkpointing
138class Data_format:
139    """Generic interface to data formats
140    """
141
142       
143    def __init__(self, domain, extension, mode = 'w'):
144        assert mode in ['r', 'w', 'a'], '''Mode %s must be either:''' %mode +\
145                                        '''   'w' (write)'''+\
146                                        '''   'r' (read)''' +\
147                                        '''   'a' (append)'''   
148
149        #Create filename
150        self.filename = create_filename(domain.get_name(), extension, len(domain))
151        self.timestep = 0
152        self.number_of_volumes = len(domain)
153        self.domain = domain
154   
155   
156    #FIXME: Should we have a general set_precision function?
157   
158
159
160#Class for storing output to e.g. visualisation
161class Data_format_sww(Data_format):
162    """Interface to native NetCDF format (.sww)
163    """
164
165       
166    def __init__(self, domain, mode = 'w'):
167        from Scientific.IO.NetCDF import NetCDFFile
168        from Numeric import Int, Float, Float32
169
170        self.precision = Float32 #Use single precision
171       
172        Data_format.__init__(self, domain, 'sww', mode)
173
174
175        # NetCDF file definition
176        fid = NetCDFFile(self.filename, mode)
177       
178        if mode == 'w':
179           
180            #Create new file
181            fid.institution = 'Geoscience Australia'
182            fid.description = 'Output from pyvolution suitable for plotting'
183           
184            if domain.smooth:
185                fid.smoothing = 'Yes'
186            else:   
187                fid.smoothing = 'No'
188
189            fid.order = domain.default_order   
190
191            # dimension definitions
192            fid.createDimension('number_of_volumes', self.number_of_volumes) 
193            fid.createDimension('number_of_vertices', 3)
194
195            if domain.smooth is True:
196                fid.createDimension('number_of_points', len(domain.vertexlist))
197            else:   
198                fid.createDimension('number_of_points', 3*self.number_of_volumes)
199               
200            fid.createDimension('number_of_timesteps', None) #extensible
201
202            # variable definitions
203            fid.createVariable('x', self.precision, ('number_of_points',))
204            fid.createVariable('y', self.precision, ('number_of_points',))
205            fid.createVariable('z', self.precision, ('number_of_points',))
206       
207            fid.createVariable('volumes', Int, ('number_of_volumes',
208                                                'number_of_vertices'))
209
210            fid.createVariable('time', self.precision,
211                               ('number_of_timesteps',))
212       
213            fid.createVariable('stage', self.precision,
214                               ('number_of_timesteps',
215                                'number_of_points'))
216
217            fid.createVariable('uh', self.precision,
218                               ('number_of_timesteps',
219                                'number_of_points'))
220
221            fid.createVariable('vh', self.precision,
222                               ('number_of_timesteps',
223                                'number_of_points'))
224
225        #Close
226        fid.close()               
227
228
229    def store_connectivity(self):
230        """Specialisation of store_connectivity for net CDF format
231
232        Writes x,y,z coordinates of triangles constituting
233        the bed elevation.
234        """
235
236        from Scientific.IO.NetCDF import NetCDFFile
237
238        from Numeric import concatenate
239
240        domain = self.domain
241
242        #Get NetCDF
243        fid = NetCDFFile(self.filename, 'a')  #Open existing file for append
244       
245        # Get the variables
246        x = fid.variables['x']
247        y = fid.variables['y']
248        z = fid.variables['z']
249       
250        volumes = fid.variables['volumes']
251
252        # Get X, Y and bed elevation Z
253        Q = domain.quantities['elevation']
254        X,Y,Z,V = Q.get_vertex_values(xy=True,
255                                      precision = self.precision)
256                                     
257
258
259        x[:] = X.astype(self.precision)
260        y[:] = Y.astype(self.precision)
261        z[:] = Z.astype(self.precision)
262
263        volumes[:] = V
264       
265        #Close
266        fid.close()
267
268
269    def store_timestep(self, names):
270        """Store time and named quantities to file
271        """
272        from Scientific.IO.NetCDF import NetCDFFile
273        import types
274        from time import sleep
275
276               
277        #Get NetCDF
278        retries = 0
279        file_open = False
280        while not file_open and retries < 10:
281            try:       
282                fid = NetCDFFile(self.filename, 'a') #Open existing file
283            except IOError:
284                #This could happen if someone was reading the file.
285                #In that case, wait a while and try again
286                msg = 'Warning (store_timestep): File %s could not be opened'\
287                %self.filename
288                msg += ' - trying again'
289                print msg
290                retries += 1
291                sleep(1)
292            else:
293               file_open = True
294               
295        if not file_open:
296            msg = 'File %s could not be opened for append' %self.filename
297            raise msg
298               
299           
300        domain = self.domain
301       
302        # Get the variables
303        time = fid.variables['time']
304        stage = fid.variables['stage']
305        uh = fid.variables['uh']
306        vh = fid.variables['vh']                       
307        i = len(time)
308
309        #Store time
310        time[i] = self.domain.time
311
312
313        if type(names) not in [types.ListType, types.TupleType]:
314            names = [names]
315
316        for name in names:   
317            # Get quantity
318            Q = domain.quantities[name]
319            A,V = Q.get_vertex_values(xy=False,
320                                      precision = self.precision)
321
322            #FIXME: Make this more general and rethink naming
323            if name == 'level':
324                stage[i,:] = A.astype(self.precision)
325            elif name == 'xmomentum':
326                uh[i,:] = A.astype(self.precision)
327            elif name == 'ymomentum':
328                vh[i,:] = A.astype(self.precision)                       
329           
330
331        #Flush and close
332        fid.sync()
333        fid.close()
334
335
336#Class for handling checkpoints data
337class Data_format_cpt(Data_format):
338    """Interface to native NetCDF format (.cpt)
339    """
340
341       
342    def __init__(self, domain, mode = 'w'):
343        from Scientific.IO.NetCDF import NetCDFFile
344        from Numeric import Int, Float, Float
345
346        self.precision = Float #Use full precision
347       
348        Data_format.__init__(self, domain, 'sww', mode)
349
350
351        # NetCDF file definition
352        fid = NetCDFFile(self.filename, mode)
353       
354        if mode == 'w':
355            #Create new file
356            fid.institution = 'Geoscience Australia'
357            fid.description = 'Checkpoint data'
358            #fid.smooth = domain.smooth               
359            fid.order = domain.default_order   
360
361            # dimension definitions
362            fid.createDimension('number_of_volumes', self.number_of_volumes) 
363            fid.createDimension('number_of_vertices', 3)
364
365            #Store info at all vertices (no smoothing)
366            fid.createDimension('number_of_points', 3*self.number_of_volumes)
367            fid.createDimension('number_of_timesteps', None) #extensible
368
369            # variable definitions
370
371            #Mesh
372            fid.createVariable('x', self.precision, ('number_of_points',))
373            fid.createVariable('y', self.precision, ('number_of_points',))
374            #fid.createVariable('z', self.precision, ('number_of_points',))
375       
376            fid.createVariable('volumes', Int, ('number_of_volumes',
377                                                'number_of_vertices'))
378
379            fid.createVariable('time', self.precision,
380                               ('number_of_timesteps',))
381
382            #Allocate space for all quantities
383            for name in domain.quantities.keys():
384                fid.createVariable(name, self.precision,
385                                   ('number_of_timesteps',
386                                    'number_of_points'))
387
388        #Close
389        fid.close()               
390
391
392    def store_checkpoint(self):
393        """
394        Write x,y coordinates of triangles.
395        Write connectivity (
396        constituting
397        the bed elevation.
398        """
399
400        from Scientific.IO.NetCDF import NetCDFFile
401
402        from Numeric import concatenate
403
404        domain = self.domain
405
406        #Get NetCDF
407        fid = NetCDFFile(self.filename, 'a')  #Open existing file for append
408       
409        # Get the variables
410        x = fid.variables['x']
411        y = fid.variables['y']
412       
413        volumes = fid.variables['volumes']
414
415        # Get X, Y and bed elevation Z
416        Q = domain.quantities['elevation']
417        X,Y,Z,V = Q.get_vertex_values(xy=True,
418                                      precision = self.precision)
419                                     
420
421
422        x[:] = X.astype(self.precision)
423        y[:] = Y.astype(self.precision)
424        z[:] = Z.astype(self.precision)
425
426        volumes[:] = V
427       
428        #Close
429        fid.close()
430
431
432    def store_timestep(self, name):
433        """Store time and named quantity to file
434        """
435        from Scientific.IO.NetCDF import NetCDFFile
436        from time import sleep
437               
438        #Get NetCDF
439        retries = 0
440        file_open = False
441        while not file_open and retries < 10:
442            try:       
443                fid = NetCDFFile(self.filename, 'a') #Open existing file
444            except IOError:
445                #This could happen if someone was reading the file.
446                #In that case, wait a while and try again
447                msg = 'Warning (store_timestep): File %s could not be opened'\
448                %self.filename
449                msg += ' - trying again'
450                print msg
451                retries += 1
452                sleep(1)
453            else:
454               file_open = True
455               
456        if not file_open:
457            msg = 'File %s could not be opened for append' %self.filename
458            raise msg
459               
460           
461        domain = self.domain
462       
463        # Get the variables
464        time = fid.variables['time']
465        stage = fid.variables['stage']       
466        i = len(time)
467
468        #Store stage
469        time[i] = self.domain.time
470
471        # Get quantity
472        Q = domain.quantities[name]
473        A,V = Q.get_vertex_values(xy=False,
474                                  precision = self.precision)
475       
476        stage[i,:] = A.astype(self.precision)
477
478        #Flush and close
479        fid.sync()
480        fid.close()
481
482
483
484
485
486
487
488#Function for storing xya output
489#FIXME Not done yet for this version
490class Data_format_xya(Data_format):
491    """Generic interface to data formats
492    """
493
494
495    def __init__(self, domain, mode = 'w'):
496        from Scientific.IO.NetCDF import NetCDFFile
497        from Numeric import Int, Float, Float32
498
499        self.precision = Float32 #Use single precision
500       
501        Data_format.__init__(self, domain, 'xya', mode)
502       
503   
504   
505
506    def store_all(self):
507        """Specialisation of store all for xya format
508
509        Writes x,y,z coordinates of triangles constituting
510        the bed elevation.
511        """
512
513        from Numeric import concatenate
514
515        domain = self.domain
516       
517        fd = open(self.filename, 'w')
518
519
520        if domain.smooth is True:
521            number_of_points =  len(domain.vertexlist)
522        else:   
523            number_of_points = 3*self.number_of_volumes
524
525        numVertAttrib = 3 #Three attributes is what is assumed by the xya format
526
527        fd.write(str(number_of_points) + " " + str(numVertAttrib) +\
528                 " # <vertex #> <x> <y> [attributes]" + "\n")
529
530
531        # Get X, Y, bed elevation and friction (index=0,1)
532        X,Y,A,V = domain.get_vertex_values(xy=True, value_array='field_values',
533                                           indices = (0,1), precision = self.precision)
534
535        bed_eles = A[:,0]
536        fricts = A[:,1]
537
538        # Get stage (index=0)
539        B,V = domain.get_vertex_values(xy=False, value_array='conserved_quantities',
540                                       indices = (0,), precision = self.precision)
541
542        stages = B[:,0]
543
544        #<vertex #> <x> <y> [attributes]
545        for x, y, bed_ele, stage, frict in map(None, X, Y, bed_eles,
546                                               stages, fricts):
547
548            s = '%.6f %.6f %.6f %.6f %.6f\n' %(x, y, bed_ele, stage, frict)
549            fd.write(s)
550           
551        #close
552        fd.close()
553
554
555    def store_timestep(self, t, V0, V1, V2):
556        """Store time, water heights (and momentums) to file
557        """
558        pass
559   
560
561#Auxiliary
562def write_obj(filename,x,y,z):
563    """Store x,y,z vectors into filename (obj format)
564       Vectors are assumed to have dimension (M,3) where
565       M corresponds to the number elements.
566       triangles are assumed to be disconnected
567
568       The three numbers in each vector correspond to three vertices,
569
570       e.g. the x coordinate of vertex 1 of element i is in x[i,1]
571       
572    """
573    #print 'Writing obj to %s' % filename
574
575    import os.path
576
577    root, ext = os.path.splitext(filename)
578    if ext == '.obj':
579        FN = filename
580    else:
581        FN = filename + '.obj'
582       
583
584    outfile = open(FN, 'wb')
585    outfile.write("# Triangulation as an obj file\n")
586
587    M, N = x.shape
588    assert N==3  #Assuming three vertices per element
589
590    for i in range(M):
591        for j in range(N):
592            outfile.write("v %f %f %f\n" % (x[i,j],y[i,j],z[i,j]))
593
594    for i in range(M):
595        base = i*N
596        outfile.write("f %d %d %d\n" % (base+1,base+2,base+3))
597
598    outfile.close()
599
600
601
602#Conversion routines
603def sww2obj(basefilename, size):
604    """Convert netcdf based data output to obj
605    """
606    from Scientific.IO.NetCDF import NetCDFFile       
607
608    from Numeric import Float, zeros
609
610    #Get NetCDF
611    FN = create_filename(basefilename, 'sww', size)
612    print 'Reading from ', FN
613    fid = NetCDFFile(FN, 'r')  #Open existing file for read
614
615       
616    # Get the variables
617    x = fid.variables['x']
618    y = fid.variables['y']
619    z = fid.variables['z']       
620    time = fid.variables['time']
621    stage = fid.variables['stage']
622       
623    M = size  #Number of lines
624    xx = zeros((M,3), Float)
625    yy = zeros((M,3), Float)
626    zz = zeros((M,3), Float)
627
628    for i in range(M):
629        for j in range(3):
630            xx[i,j] = x[i+j*M]
631            yy[i,j] = y[i+j*M]
632            zz[i,j] = z[i+j*M]
633       
634    #Write obj for bathymetry
635    FN = create_filename(basefilename, 'obj', size)       
636    write_obj(FN,xx,yy,zz)
637
638
639    #Now read all the data with variable information, combine with
640    #x,y info and store as obj
641
642    for k in range(len(time)):
643        t = time[k]
644        print 'Processing timestep %f' %t
645
646        for i in range(M):
647            for j in range(3):
648                zz[i,j] = stage[k,i+j*M]
649
650
651        #Write obj for variable data
652        #FN = create_filename(basefilename, 'obj', size, time=t)
653        FN = create_filename(basefilename[:5], 'obj', size, time=t)   
654        write_obj(FN,xx,yy,zz)
655
656
657def dat2obj(basefilename):
658    """Convert line based data output to obj
659    """
660
661    import glob, os
662    from config import data_dir
663   
664   
665    #Get bathymetry and x,y's
666    lines = open(data_dir+os.sep+basefilename+'_geometry.dat', 'r').readlines()
667   
668    from Numeric import zeros, Float
669
670    M = len(lines)  #Number of lines
671    x = zeros((M,3), Float)
672    y = zeros((M,3), Float)
673    z = zeros((M,3), Float)
674
675    ##i = 0
676    for i, line in enumerate(lines):
677        tokens = line.split()
678        values = map(float,tokens)       
679
680        for j in range(3):
681            x[i,j] = values[j*3]
682            y[i,j] = values[j*3+1]
683            z[i,j] = values[j*3+2]       
684       
685        ##i += 1
686
687
688    #Write obj for bathymetry
689    write_obj(data_dir+os.sep+basefilename+'_geometry',x,y,z)
690
691
692    #Now read all the data files with variable information, combine with
693    #x,y info
694    #and store as obj
695
696    files = glob.glob(data_dir+os.sep+basefilename+'*.dat')
697
698    for filename in files:
699        print 'Processing %s' % filename
700
701        lines = open(data_dir+os.sep+filename,'r').readlines()
702        assert len(lines) == M
703        root, ext = os.path.splitext(filename)
704
705        #Get time from filename
706        i0 = filename.find('_time=')
707        if i0 == -1:
708            #Skip bathymetry file
709            continue
710       
711        i0 += 6  #Position where time starts           
712        i1 = filename.find('.dat')         
713
714        if i1 > i0:
715            t = float(filename[i0:i1])
716        else:
717            raise 'Hmmmm'
718               
719       
720       
721        ##i = 0
722        for i, line in enumerate(lines):
723            tokens = line.split()
724            values = map(float,tokens)       
725
726            for j in range(3):
727                z[i,j] = values[j]       
728       
729            ##i += 1
730
731        #Write obj for variable data
732        write_obj(data_dir+os.sep+basefilename+'_time=%.4f' %t,x,y,z)
733
734
735def filter_netcdf(filename1, filename2, first=0, last=None, step = 1):
736    """Read netcdf filename1, pick timesteps first:step:last and save to
737    nettcdf file filename2
738    """
739    from Scientific.IO.NetCDF import NetCDFFile       
740   
741    #Get NetCDF
742    infile = NetCDFFile(filename1, 'r')  #Open existing file for read
743    outfile = NetCDFFile(filename2, 'w')  #Open new file
744   
745
746    #Copy dimensions
747    for d in infile.dimensions:
748        outfile.createDimension(d, infile.dimensions[d])
749
750    for name in infile.variables:
751        var = infile.variables[name]
752        outfile.createVariable(name, var.typecode(), var.dimensions)
753       
754   
755    #Copy the static variables
756    for name in infile.variables:
757        if name == 'time' or name == 'stage':
758            pass
759        else:
760            #Copy
761            outfile.variables[name][:] = infile.variables[name][:]           
762
763    #Copy selected timesteps
764    time = infile.variables['time']
765    stage = infile.variables['stage']
766
767    newtime = outfile.variables['time']
768    newstage = outfile.variables['stage']   
769
770    if last is None:
771        last = len(time)
772
773    selection = range(first, last, step)
774    for i, j in enumerate(selection):
775        print 'Copying timestep %d of %d (%f)' %(j, last-first, time[j])
776        newtime[i] = time[j]
777        newstage[i,:] = stage[j,:]       
778       
779    #Close
780    infile.close()
781    outfile.close()
782
783
784#Get data objects
785def get_dataobject(domain, mode='w'):
786    """Return instance of class of given format using filename
787    """
788
789    cls = eval('Data_format_%s' %domain.format)
790    return cls(domain, mode)
791
792
793
794
795#OBSOLETE STUFF
796#Native checkpoint format.
797#Information needed to recreate a state is preserved
798#FIXME: Rethink and maybe use netcdf format
799def cpt_variable_writer(filename, t, v0, v1, v2):
800    """Store all conserved quantities to file
801    """
802
803    M, N = v0.shape
804       
805    FN = create_filename(filename, 'cpt', M, t)
806    #print 'Writing to %s' %FN
807   
808    fid = open(FN, 'w')
809    for i in range(M):
810        for j in range(N):       
811            fid.write('%.16e ' %v0[i,j])           
812        for j in range(N):
813            fid.write('%.16e ' %v1[i,j])                       
814        for j in range(N):
815            fid.write('%.16e ' %v2[i,j])                       
816           
817        fid.write('\n')
818    fid.close()
819
820
821def cpt_variable_reader(filename, t, v0, v1, v2):
822    """Store all conserved quantities to file
823    """
824
825    M, N = v0.shape
826       
827    FN = create_filename(filename, 'cpt', M, t)
828    #print 'Reading from %s' %FN
829   
830    fid = open(FN)
831
832
833    for i in range(M):
834        values = fid.readline().split() #Get one line
835
836        for j in range(N):
837            v0[i,j] = float(values[j])
838            v1[i,j] = float(values[3+j])
839            v2[i,j] = float(values[6+j])
840   
841    fid.close()
842
843def cpt_constant_writer(filename, X0, X1, X2, v0, v1, v2):
844    """Writes x,y,z,z,z coordinates of triangles constituting the bed
845    elevation.
846    Not in use pt
847    """
848
849    M, N = v0.shape
850
851    print X0
852    import sys; sys.exit() 
853    FN = create_filename(filename, 'cpt', M)
854    print 'Writing to %s' %FN
855   
856    fid = open(FN, 'w')
857    for i in range(M):
858        for j in range(2):       
859            fid.write('%.16e ' %X0[i,j])   #x, y
860        for j in range(N):                   
861            fid.write('%.16e ' %v0[i,j])       #z,z,z,
862           
863        for j in range(2):                   
864            fid.write('%.16e ' %X1[i,j])   #x, y
865        for j in range(N):                               
866            fid.write('%.16e ' %v1[i,j]) 
867       
868        for j in range(2):                               
869            fid.write('%.16e ' %X2[i,j])   #x, y
870        for j in range(N):                                           
871            fid.write('%.16e ' %v2[i,j]) 
872           
873        fid.write('\n')
874    fid.close()
875
876
877
878#Function for storing out to e.g. visualisation
879#FIXME: Do we want this?
880#FIXME: Not done yet for this version
881def dat_constant_writer(filename, X0, X1, X2, v0, v1, v2):
882    """Writes x,y,z coordinates of triangles constituting the bed elevation.
883    """
884
885    M, N = v0.shape
886       
887    FN = create_filename(filename, 'dat', M)
888    #print 'Writing to %s' %FN
889   
890    fid = open(FN, 'w')
891    for i in range(M):
892        for j in range(2):       
893            fid.write('%f ' %X0[i,j])   #x, y
894        fid.write('%f ' %v0[i,0])       #z
895           
896        for j in range(2):                   
897            fid.write('%f ' %X1[i,j])   #x, y
898        fid.write('%f ' %v1[i,0])       #z
899       
900        for j in range(2):                               
901            fid.write('%f ' %X2[i,j])   #x, y
902        fid.write('%f ' %v2[i,0])       #z
903           
904        fid.write('\n')
905    fid.close()
906
907
908
909def dat_variable_writer(filename, t, v0, v1, v2):
910    """Store water height to file
911    """
912
913    M, N = v0.shape
914       
915    FN = create_filename(filename, 'dat', M, t)
916    #print 'Writing to %s' %FN
917   
918    fid = open(FN, 'w')
919    for i in range(M):
920        fid.write('%.4f ' %v0[i,0])
921        fid.write('%.4f ' %v1[i,0])
922        fid.write('%.4f ' %v2[i,0])       
923
924        fid.write('\n')
925    fid.close()
926
Note: See TracBrowser for help on using the repository browser.