Added slice plot to the resolution matrices.
authorAndy Kass <kass.andrew@gmail.com>
Mon, 4 Sep 2017 07:23:02 +0000 (09:23 +0200)
committerAndy Kass <kass.andrew@gmail.com>
Mon, 4 Sep 2017 07:23:02 +0000 (09:23 +0200)
EAGE_Abstract/Lemma/psf/psf-4-modified.py

index a20bf3c..2e6dc05 100644 (file)
@@ -14,6 +14,7 @@ PS = []
 K1 = []
 R = []
 
+svs = 10
 for ii in range(1,5):
     kf = open(sys.argv[ii])
     ifaces = np.array (kf.readline().split(),dtype=np.float)
@@ -31,7 +32,7 @@ for ii in range(1,5):
     GTGI = np.linalg.inv(np.dot(KA.T,KA))
     R.append(np.dot(GTGI,np.dot(KA.T,KA)))
     u,s,v = np.linalg.svd(KK.T)
-    PS.append(np.dot(v[:,0:16],np.transpose(v[:,0:16])))
+    PS.append(np.dot(v[:,0:svs],np.transpose(v[:,0:svs])))
     kf.close()
 
 K1T = np.concatenate(K1,axis=1)
@@ -39,7 +40,7 @@ KAT = np.abs(K1T)
 KKT = np.dot(KAT,KAT.T)
 R.append(np.dot(np.linalg.inv(KKT),np.dot(KKT,KAT)))
 u, s, v = np.linalg.svd(KKT.T)
-PS.append(np.dot(v[:,0:16],np.transpose(v[:,0:16])))
+PS.append(np.dot(v[:,0:svs],np.transpose(v[:,0:svs])))
 
 
 #Plot SVD resolution
@@ -57,10 +58,10 @@ axcb = fig.add_axes( [0.900, 0.25, .025, 0.583] )
 im = []
 ttls = ['Coincident', 'X Comp', 'Y Comp', 'Z Comp']
 for ii in range(0,4):
-    im.append(ax[ii].pcolormesh(X,Y, PS[ii], cmap='viridis'))#, vmin=0,
-    #im.append(ax[ii].pcolormesh(X,Y, np.sqrt(PS[ii]), cmap='viridis'))#, vmin=0,
+    im.append(ax[ii].contourf(X[0:-1,0:-1],Y[0:-1,0:-1], PS[ii]))#, vmin=0,
+    #im.append(ax[ii].pcolormesh(X,Y, PS[ii], cmap='viridis'))#, vmin=0,
         #vmax=3100))
-    im[ii].set_edgecolor('face')
+    #im[ii].set_edgecolor('face')
     ax[ii].set_ylim( ifaces[-1], ifaces[0] )
     ax[ii].set_xlim( ifaces[0], ifaces[-1] )
     ax[ii].set_xlabel('Depth [m]')
@@ -103,6 +104,16 @@ ax2.set_xlabel('Depth [m]')
 plt.show()
 plt.close()
 
+#Grab a slice from the resolution matrix
+r = 10
+fig = plt.figure()
+ax = plt.plot(X[0:-1,r],PS[0][:,r])
+plt.show()
+plt.close()
+
+
+
+
 #Plot explicit resolution 
 print(R[0].shape)
 fig = plt.figure ( figsize=(pc2in(80),pc2in(24)) )
@@ -135,7 +146,7 @@ for a in ax[1:]:
 cb = plt.colorbar(im[0], axcb )
 cb.set_label(r'$\left| \mathcal{V}_0 \right|$~(nV)')
 
-plt.savefig('RES_R.pdf',bbox_inches='tight')
-plt.show()
+#plt.savefig('RES_R.pdf',bbox_inches='tight')
+#plt.show()
 plt.close()