Changeset 9057

Timestamp:

Mar 3, 2014, 10:05:11 AM (11 years ago)

Author:

davies

Message:

Fixing towradgi test following Petar's suggestion

Location:

trunk/anuga_core/source/anuga_validation_tests/case_studies/towradgi

Files:

• Compare_results_with_fieldObs.py (modified) (5 diffs)
• Validation (added)
• Validation/historic_1998_flood_levels_towradgi_ck.csv (added)
• Validation/pioneer_timeseries.txt (added)

trunk/anuga_core/source/anuga_validation_tests/case_studies/towradgi/Compare_results_with_fieldObs.py

@@ -8 +8 @@
 from anuga.utilities import plot_utils as util
 
+
 swwdir='MODEL_OUTPUTS/'
-#swwdir='MODEL_bridges_n0p04/'
+swwname='Towradgi_historic_flood.sww' #'Aug98_hort_discontinous.sww' #'Towradgi_historic_flood.sww'
 
-p=util.get_output(swwdir+'Towradgi_historic_flood.sww')
+p=util.get_output(swwdir+swwname)
 p2=util.get_centroids(p,velocity_extrapolation=True)
 floodLevels=scipy.genfromtxt('Validation/historic_1998_flood_levels_towradgi_ck.csv', delimiter=',',skip_header=1)
+pioneerLevel=scipy.genfromtxt('Validation/pioneer_timeseries.txt', 
+                              skip_header=1)
 
 # Extract modelled peak at the coordinates in floodLevels
@@ -27 +30 @@
     modelled_ind[i]=myInd
     modelled_level[i]=p2.stage[:,myInd].max()
+
+    if(i==0):
+        pyplot.clf()
+        pyplot.plot(p2.time,p2.stage[:, myInd],color='red')
+        pyplot.plot((pioneerLevel[:,0]-540.)*60., pioneerLevel[:,2],color='black')
+        pyplot.xlabel('Time')
+        pyplot.ylabel('Stage (m)')
+        pyplot.title('Pioneer Bridge Stage')
+        pyplot.savefig('Pioneer_Bridge_Stage.png')
+    
 
 pyplot.clf()
@@ -53 +66 @@
     tif_outdir='OUTPUT_TIFS'
     CellSize=5.0
-    util.Make_Geotif(swwdir+'Towradgi_historic_flood.sww',
+    print 'Making tifs'
+    util.Make_Geotif(swwdir+swwname,
                       ['depth','velocity','depthIntegratedVelocity','elevation'],'max',
                       CellSize=CellSize,EPSG_CODE=32756,output_dir=tif_outdir)
-
+    print 'Made tifs'
     # Plot depth raster with discrepency between model and data
     depthFile=tif_outdir+'/Towradgi_historic_flood_depth_max_max.tif'
@@ -66 +80 @@
     pyplot.figure(figsize=(12,6))
     pyplot.plot([X.min(),X.max()],[Y.min(),Y.max()],' ')
-    pyplot.imshow(myDepth,extent=[X.min(),X.max(),Y.min(),Y.max()],origin='lower',cmap=pyplot.get_cmap('Greys'))
+    print 'Imshow'
+    #import pdb
+    #pdb.set_trace()
+    pyplot.imshow(scipy.flipud(myDepth),extent=[X.min(),X.max(),Y.min(),Y.max()],origin='lower',cmap=pyplot.get_cmap('Greys'))
     pyplot.gca().set_aspect('equal')
     pyplot.colorbar(orientation='horizontal').set_label('Peak Depth in model (m)')
-    pyplot.scatter(floodLevels[:,0], floodLevels[:,1], c=floodLevels[:,3]-modelled_level,s=20,cmap=pyplot.get_cmap('spectral'))
+    er1=floodLevels[:,3]-modelled_level
+    #er1=er1*(er1<1.0) + 1.0*(er1>=1.0)
+    #er1=er1*(er1> -1.0) - 1.0*(er1<=-1.0)
+    pyplot.scatter(floodLevels[:,0], floodLevels[:,1], c=er1,s=20,cmap=pyplot.get_cmap('spectral'))
     pyplot.colorbar().set_label(label='Field observation - Modelled Peak Depth (m)')
     pyplot.xlim([p.x.min()+p.xllcorner,p.x.max()+p.xllcorner])
@@ -75 +95 @@
     pyplot.savefig('Spatial_Depth_and_Error.png')
 except:
-    'Cannot make GIS plot -- perhaps GDAL etc are not installed?'
+    print 'Cannot make GIS plot -- perhaps GDAL etc are not installed?'