Author: bugman Date: Tue Oct 14 13:41:57 2014 New Revision: 26276 URL: http://svn.gna.org/viewcvs/relax?rev=26276&view=rev Log: Merged revisions 26256-26258 via svnmerge from svn+ssh://bugman@xxxxxxxxxxx/svn/relax/trunk ........ r26256 | bugman | 2014-10-13 16:30:13 +0200 (Mon, 13 Oct 2014) | 6 lines The devel_scripts/python_multiversion_test_suite.py script now runs relax with the --time flag. This is for quicker identification of failure points. It will also force the sys.stdout buffer to be flushed more often on Python 2.5 so that it does not appear as if the tests have frozen. ........ r26257 | tlinnet | 2014-10-13 17:18:50 +0200 (Mon, 13 Oct 2014) | 1 line Added check to systemtest Relax_disp.test_cpmg_synthetic_dx_map_points() for the creation of a matplotlib surface command plot file. ........ r26258 | tlinnet | 2014-10-13 17:18:56 +0200 (Mon, 13 Oct 2014) | 10 lines Added the write out of a matplotlib command file, to plot surfaces of a dx map. It uses the minimum chi2 value in the map space, to define surface defitions. It creates a X,Y; X,Z; Y,Z map, where the values in the missing dimension has been cut at the minimum chi2 value. For each map, it creates a projected 3d map of the parameters and the chi2 value, and a heat map for the contours. It also scatters the minimum chi2 value, the 4 smallest chi2 values, and maps any points in the point file, to a scatter point. Mapping the points from file to map points, is done by finding the shortest Euclidean distance in the space from the points to any map points. ........ Modified: branches/frame_order_cleanup/ (props changed) branches/frame_order_cleanup/devel_scripts/python_multiversion_test_suite.py branches/frame_order_cleanup/pipe_control/opendx.py branches/frame_order_cleanup/test_suite/system_tests/relax_disp.py Propchange: branches/frame_order_cleanup/ ------------------------------------------------------------------------------ --- svnmerge-integrated (original) +++ svnmerge-integrated Tue Oct 14 13:41:57 2014 @@ -1 +1 @@ -/trunk:1-26205,26208-26255 +/trunk:1-26205,26208-26258 Modified: branches/frame_order_cleanup/devel_scripts/python_multiversion_test_suite.py URL: http://svn.gna.org/viewcvs/relax/branches/frame_order_cleanup/devel_scripts/python_multiversion_test_suite.py?rev=26276&r1=26275&r2=26276&view=diff ============================================================================== --- branches/frame_order_cleanup/devel_scripts/python_multiversion_test_suite.py (original) +++ branches/frame_order_cleanup/devel_scripts/python_multiversion_test_suite.py Tue Oct 14 13:41:57 2014 @@ -85,4 +85,4 @@ execute_sh('~/bin/python%s relax -i' % (version), log=LOG, err=LOG) # Run the test suite. - execute_sh('~/bin/python%s relax -x' % (version), log=LOG, err=LOG) + execute_sh('~/bin/python%s relax -x --time' % (version), log=LOG, err=LOG) Modified: branches/frame_order_cleanup/pipe_control/opendx.py URL: http://svn.gna.org/viewcvs/relax/branches/frame_order_cleanup/pipe_control/opendx.py?rev=26276&r1=26275&r2=26276&view=diff ============================================================================== --- branches/frame_order_cleanup/pipe_control/opendx.py (original) +++ branches/frame_order_cleanup/pipe_control/opendx.py Tue Oct 14 13:41:57 2014 @@ -174,6 +174,9 @@ if create_par_file: self.create_par_chi2(file_prefix=self.file_prefix, par_chi2_vals=self.par_chi2_vals) + # Generate the matplotlib script file, to plot surfaces + self.matplotlib_surface_plot() + ## Generate the file with parameters and associated chi2 value for the points send to dx. if self.num_points >= 1 and create_par_file: # Calculate the parameter and associated chi2 values for the points. @@ -449,3 +452,235 @@ string = string + " " + repr(val) val = val + loc_inc self.tick_locations.append("{" + string + " }") + + + def matplotlib_surface_plot(self): + """Function to write matplotlib script to plot surfaces of parameters.""" + + # Add ".par" to file prefix + mapfile_name = '"%s.par"' % self.file_prefix + + # If point file_file is different from None + if self.point_file != None: + pointfile_name = '"%s.par"' % self.point_file + else: + pointfile_name = "None" + + # Open the file. + plot_file = open_write_file(file_name=self.file_prefix+'.py', dir=self.dir, force=True) + + matplotlib_file = [ + 'import numpy as np'+"\n", + 'import scipy.interpolate'+"\n", + 'from numpy.ma import masked_where'+"\n", + ''+"\n", + 'from mpl_toolkits.mplot3d import axes3d'+"\n", + 'import matplotlib.pyplot as plt'+"\n", + 'from matplotlib import cm'+"\n", + ''+"\n", + '# Open file and get header.'+"\n", + 'mapfile_name = %s'%mapfile_name+"\n", + 'pointfile_name = %s'%pointfile_name+"\n", + ''+"\n", + 'mapfile = open(mapfile_name, "r")'+"\n", + 'lines = mapfile.readlines()'+"\n", + 'mapfile.close()'+"\n", + 'header = lines[0].split()[1:]'+"\n", + ''+"\n", + '# Prepare the dtype for reading file.'+"\n", + 'dtype_str = "i8,f8,f8,f8,f8,i8,f8,f8,f8,f8"'+"\n", + ''+"\n", + 'print("Fileheader is: %s"%header)'+"\n", + 'print("Value types are: %s"%dtype_str)'+"\n", + ''+"\n", + '# Load the data.'+"\n", + 'data = np.genfromtxt(fname=mapfile_name, dtype=dtype_str, names=header)'+"\n", + ''+"\n", + '# Load the point data'+"\n", + 'if pointfile_name:'+"\n", + ' # Load the point data.'+"\n", + ' data_p = np.genfromtxt(fname=pointfile_name, dtype=dtype_str, names=header)'+"\n", + ' '+"\n", + ''+"\n", + '# Define where to cut the data, as the minimum.'+"\n", + 'header_min = header[6:10]'+"\n", + ''+"\n", + '# Define to cut at min map point.'+"\n", + 'map_min_par0 = data[header_min[0]][0]'+"\n", + 'map_min_par1 = data[header_min[1]][0]'+"\n", + 'map_min_par2 = data[header_min[2]][0]'+"\n", + 'map_min_chi2 = data[header_min[3]][0]'+"\n", + ''+"\n", + '# Now get the headers for the data.'+"\n", + 'header_val = header[1:5]'+"\n", + ''+"\n", + '# Define which 2D maps to create, as a list of 2 parameters, and at which third parameter to cut the values.'+"\n", + 'maps_xy = [header_val[0], header_val[1], header_val[2], map_min_par2]'+"\n", + 'maps_xz = [header_val[0], header_val[2], header_val[1], map_min_par1]'+"\n", + 'maps_yz = [header_val[1], header_val[2], header_val[0], map_min_par0]'+"\n", + ''+"\n", + 'maps = [maps_xy, maps_xz, maps_yz]'+"\n", + ''+"\n", + '# Nr of columns is number of maps.'+"\n", + 'nr_cols = 1'+"\n", + '# Nr of rows, is 2, for 3d projection and imshow'+"\n", + 'nr_rows = 2'+"\n", + ''+"\n", + '# Loop over the maps:'+"\n", + 'for x_par, y_par, z_par, z_cut in maps:'+"\n", + ' # Define figure'+"\n", + ' fig = plt.figure()'+"\n", + ''+"\n", + ' # Define c_par'+"\n", + ' c_par = header_val[3]'+"\n", + ''+"\n", + ' # Now get the values for the map.'+"\n", + ' map_x = data[x_par]'+"\n", + ' map_y = data[y_par]'+"\n", + ' map_z = data[z_par]'+"\n", + ' map_c = data[c_par]'+"\n", + ''+"\n", + ' # Now define which map to create.'+"\n", + ' mask_xy = masked_where(map_z == z_cut, map_z)'+"\n", + ' map_mask_x = map_x[mask_xy.mask]'+"\n", + ' map_mask_y = map_y[mask_xy.mask]'+"\n", + ' map_mask_c = map_c[mask_xy.mask]'+"\n", + ''+"\n", + ' # Define min and max values.'+"\n", + ' map_mask_x_min = map_mask_x.min()'+"\n", + ' map_mask_x_max = map_mask_x.max()'+"\n", + ' map_mask_y_min = map_mask_y.min()'+"\n", + ' map_mask_y_max = map_mask_y.max()'+"\n", + ' map_mask_c_min = map_mask_c.min()'+"\n", + ' map_mask_c_max = map_mask_c.max()'+"\n", + ''+"\n", + ' # Set up a regular grid of interpolation points'+"\n", + ' int_points = 300'+"\n", + ' xi, yi = np.linspace(map_mask_x_min, map_mask_x_max, int_points), np.linspace(map_mask_y_min, map_mask_y_max, int_points)'+"\n", + ' xi, yi = np.meshgrid(xi, yi)'+"\n", + ''+"\n", + ' # Interpolate to create grid'+"\n", + ' ci = scipy.interpolate.griddata((map_mask_x, map_mask_y), map_mask_c, (xi, yi), method="linear")'+"\n", + ''+"\n", + ' # Set which x, y, z to plot'+"\n", + ' x_p = xi'+"\n", + ' y_p = yi'+"\n", + ' c_p = ci'+"\n", + ''+"\n", + ' # Cut map at a certain height.'+"\n", + ' # First get index os largest values'+"\n", + ' #out_val = 5*map_mask_c_min'+"\n", + ' out_val = map_mask_c_max'+"\n", + ' ci_mask = masked_where(ci >= out_val, ci)'+"\n", + ''+"\n", + ' # Replace with 0.0'+"\n", + ' ci[ci_mask.mask] = 0.0'+"\n", + ' # Find new max'+"\n", + ' new_max = np.max(ci)'+"\n", + ''+"\n", + ' # Insert values in array.'+"\n", + ' ci[ci_mask.mask] = new_max'+"\n", + ''+"\n", + ' # Create figure and plot'+"\n", + ' ax = fig.add_subplot(nr_rows, nr_cols, 1, projection="3d")'+"\n", + ' ax.plot_surface(x_p, y_p, c_p, rstride=8, cstride=8, alpha=0.3)'+"\n", + ''+"\n", + ' # Possible add scatter points for map.'+"\n", + ' #ax.scatter(map_x, map_y, map_c, c="b", marker="o", s=5)'+"\n", + ''+"\n", + ' # One could also make the mesh just from the values, but this require much memory.'+"\n", + ' ##ax.scatter(x_p, y_p, c_p, c="y", marker="o", s=5)'+"\n", + ''+"\n", + ' # Add contour levels on sides.'+"\n", + ' ax.contour(x_p, y_p, c_p, zdir="z", offset=0, cmap=cm.coolwarm)'+"\n", + ' ax.contour(x_p, y_p, c_p, zdir="x", offset=map_mask_x_min, cmap=cm.coolwarm)'+"\n", + ' ax.contour(x_p, y_p, c_p, zdir="y", offset=map_mask_y_min, cmap=cm.coolwarm)'+"\n", + ''+"\n", + ' # Add scatter values, for 5 lowest values.'+"\n", + ' x_par_min = x_par + "_sort"'+"\n", + ' y_par_min = y_par + "_sort"'+"\n", + ' c_par_min = c_par + "_sort"'+"\n", + ' mp_x = data[x_par_min][0:5]'+"\n", + ' mp_y = data[y_par_min][0:5]'+"\n", + ' mp_c = data[c_par_min][0:5]'+"\n", + ' ax.scatter(mp_x[0], mp_y[0], mp_c[0], c="r", marker="o", s=200)'+"\n", + ' ax.scatter(mp_x[1:], mp_y[1:], mp_c[1:], c="g", marker="o", s=100)'+"\n", + ''+"\n", + ' # Add points from file, as the closest point in map.'+"\n", + ' if pointfile_name:'+"\n", + ' if data_p[x_par].ndim == 0:'+"\n", + ' points_x = np.asarray([data_p[x_par]])'+"\n", + ' points_y = np.asarray([data_p[y_par]])'+"\n", + ' else:'+"\n", + ' points_x = data_p[x_par]'+"\n", + ' points_y = data_p[y_par]'+"\n", + ''+"\n", + ' # Normalize, by division of largest number of map'+"\n", + ' points_x_norm = points_x / map_mask_x_max'+"\n", + ' points_y_norm = points_y / map_mask_y_max'+"\n", + ' map_mask_x_norm = map_mask_x / map_mask_x_max'+"\n", + ' map_mask_y_norm = map_mask_y / map_mask_y_max'+"\n", + ''+"\n", + ' p_x = []'+"\n", + ' p_y = []'+"\n", + ' p_c = []'+"\n", + ' # Now calculate the Euclidean distance in the space, to find map point best represents the point.'+"\n", + ' for i, point_x_norm in enumerate(points_x_norm):'+"\n", + ' point_y_norm = points_y_norm[i]'+"\n", + ''+"\n", + ' # Get the distance.'+"\n", + ' dist = np.sqrt( (map_mask_x_norm - point_x_norm)**2 + (map_mask_y_norm - point_y_norm)**2)'+"\n", + ''+"\n", + ' # Return the indices of the minimum values along an axis.'+"\n", + ' min_index = np.argmin(dist)'+"\n", + ' p_x.append(map_mask_x[min_index])'+"\n", + ' p_y.append(map_mask_y[min_index])'+"\n", + ' p_c.append(map_mask_c[min_index])'+"\n", + ''+"\n", + ' # Convert to numpy array'+"\n", + ' p_x = np.asarray(p_x)'+"\n", + ' p_y = np.asarray(p_y)'+"\n", + ' p_c = np.asarray(p_c)'+"\n", + ''+"\n", + ' # Plot points'+"\n", + ' ax.scatter(p_x, p_y, p_c, c="m", marker="o", s=50)'+"\n", + ''+"\n", + ''+"\n", + ' # Set label'+"\n", + ' ax.set_xlabel("%s"%x_par)'+"\n", + ' ax.set_ylabel("%s"%y_par)'+"\n", + ' ax.set_zlabel("%s"%c_par)'+"\n", + ''+"\n", + ' # Create figure and plot'+"\n", + ' ax = fig.add_subplot(nr_rows, nr_cols, 2)'+"\n", + ' fig_imshow = ax.imshow(ci, vmin=map_mask_c_min, vmax=map_mask_c_max, origin="lower", extent=[map_mask_x_min, map_mask_x_max, map_mask_y_min, map_mask_y_max])'+"\n", + ''+"\n", + ' # Add scatter values, for 5 lowest values.'+"\n", + ' ax.scatter(mp_x[0], mp_y[0], c=mp_c[0], marker="o", s=200)'+"\n", + ' ax.scatter(mp_x[1:], mp_y[1:], c="g", marker="o", s=100)'+"\n", + ''+"\n", + ' # Also add point to this map.'+"\n", + ' if pointfile_name:'+"\n", + ' # Plot points'+"\n", + ' ax.scatter(p_x, p_y, c="m", marker="o", s=50)'+"\n", + ''+"\n", + ' # Set label'+"\n", + ' ax.set_xlabel("%s"%x_par)'+"\n", + ' ax.set_ylabel("%s"%y_par)'+"\n", + ''+"\n", + ' # Add colorbar.'+"\n", + ' fig.subplots_adjust(right=0.8)'+"\n", + ' cbar_ax = fig.add_axes([0.85, 0.15, 0.05, 0.3])'+"\n", + ' fig.colorbar(fig_imshow, cax=cbar_ax)'+"\n", + ''+"\n", + '# Show plot.'+"\n", + 'plt.show()'+"\n", + ''+"\n", + ] + + # Loop over the lines and write. + for line in matplotlib_file: + plot_file.write(line) + + # Close the file. + plot_file.close() Modified: branches/frame_order_cleanup/test_suite/system_tests/relax_disp.py URL: http://svn.gna.org/viewcvs/relax/branches/frame_order_cleanup/test_suite/system_tests/relax_disp.py?rev=26276&r1=26275&r2=26276&view=diff ============================================================================== --- branches/frame_order_cleanup/test_suite/system_tests/relax_disp.py (original) +++ branches/frame_order_cleanup/test_suite/system_tests/relax_disp.py Tue Oct 14 13:41:57 2014 @@ -2797,6 +2797,7 @@ map_net = ds.tmpdir+sep+file_name_map+".net" map_general = ds.tmpdir+sep+file_name_map+".general" map_par = get_file_path(file_name=file_name_map+".par", dir=ds.tmpdir) + map_plot = get_file_path(file_name=file_name_map+".py", dir=ds.tmpdir) point_general = ds.tmpdir+sep+file_name_point+".general" point_point = ds.tmpdir+sep+file_name_point @@ -2807,6 +2808,7 @@ self.assert_(access(map_net, F_OK)) self.assert_(access(map_general, F_OK)) self.assert_(access(map_par, F_OK)) + self.assert_(access(map_plot, F_OK)) self.assert_(access(point_general, F_OK)) self.assert_(access(point_point, F_OK)) self.assert_(access(point_par, F_OK)) @@ -2908,13 +2910,34 @@ # continue self.assertEqual(res_file[i], lines[i]) - print("\nChecking the dx point point par file.") + print("\nChecking the dx point par file.") res_file = [ '# i dw pA kex chi2 i_sort dw_sort pA_sort kex_sort chi2_sort '+"\n", '0 2.00000 0.99000 1000.00000 6185.84926 0 2.00000 0.99000 1000.00000 6185.84926 '+"\n", '1 1.92453 0.98961 1034.72206 6396.02770 1 1.92453 0.98961 1034.72206 6396.02770 '+"\n", ] file = open(point_par, 'r') + lines = file.readlines() + file.close() + for i in range(len(res_file)): + self.assertEqual(res_file[i], lines[i]) + + print("\nChecking the matplotlib surface plot file.") + res_file = [ + 'import numpy as np'+"\n", + 'import scipy.interpolate'+"\n", + 'from numpy.ma import masked_where'+"\n", + ''+"\n", + 'from mpl_toolkits.mplot3d import axes3d'+"\n", + 'import matplotlib.pyplot as plt'+"\n", + 'from matplotlib import cm'+"\n", + ''+"\n", + '# Open file and get header.'+"\n", + 'mapfile_name = "%s.par"'%file_name_map+"\n", + 'pointfile_name = "%s.par"'%file_name_point+"\n", + ''+"\n", + ] + file = open(map_plot, 'r') lines = file.readlines() file.close() for i in range(len(res_file)):