1 # Copyright (C) 2009-2012 Massimo Sandal <devicerandom@gmail.com>
2 # Pancaldi Paolo <pancaldi.paolo@gmail.com>
3 # W. Trevor King <wking@drexel.edu>
5 # This file is part of Hooke.
7 # Hooke is free software: you can redistribute it and/or modify it
8 # under the terms of the GNU Lesser General Public License as
9 # published by the Free Software Foundation, either version 3 of the
10 # License, or (at your option) any later version.
12 # Hooke is distributed in the hope that it will be useful, but WITHOUT
13 # ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14 # or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General
15 # Public License for more details.
17 # You should have received a copy of the GNU Lesser General Public
18 # License along with Hooke. If not, see
19 # <http://www.gnu.org/licenses/>.
21 from ..libhooke import WX_GOOD, ClickedPoint, config_file_path
25 wxversion.select(WX_GOOD)
26 from wx import PostEvent
34 from .. import curve as lhc
38 warnings.simplefilter('ignore',np.RankWarning)
41 class pclusterCommands(object):
47 def do_pcluster(self,args):
51 Automatically measures peaks and extracts informations for further clustering
52 (c)Paolo Pancaldi, Massimo Sandal 2009
54 blindw = str(self.convfilt_config['blindwindow'])
55 pclus_dir = "pCluster_blind"+blindw+"_"+time.strftime("%Y%m%d_%H%M")
56 self.my_work_dir = os.path.join(os.getcwd(), pclus_dir)
57 self.my_curr_dir = os.path.basename(os.getcwd())
58 os.mkdir(self.my_work_dir)
60 #--Custom persistent length
62 for arg in args.split():
63 #look for a persistent length argument.
65 pl_expression=arg.split('=')
66 pl_value=float(pl_expression[1]) #actual value
70 #configuration variables
71 min_npks = self.convfilt_config['minpeaks']
72 min_deviation = self.convfilt_config['mindeviation']
74 pclust_filename = "automeasure_"+self.my_curr_dir+"_blind"+blindw+".txt" #raw_input('Automeasure filename? ')
75 realclust_filename = "coordinate_"+self.my_curr_dir+"_blind"+blindw+".txt" #raw_input('Coordinates filename? ')
76 peackforce_filename = "peakforce_"+self.my_curr_dir+"_blind"+blindw+".txt" #raw_input('Peacks and Forces filename? ')
78 f=open(self.my_work_dir+pclust_filename,'w+')
79 f.write('Analysis started '+time.asctime()+'\n')
80 f.write('----------------------------------------\n')
81 f.write('; Contour length (nm) ; Persistence length (nm) ; Max.Force (pN) ; Slope (N/m) ; Sigma contour (nm) ; Sigma persistence (nm)\n')
84 f=open(self.my_work_dir+realclust_filename,'w+')
85 f.write('Analysis started '+time.asctime()+'\n')
86 f.write('----------------------------------------\n')
87 f.write('; Peak number ; Mean delta (nm) ; Median delta (nm) ; Mean force (pN) ; Median force (pN) ; First peak length (nm) ; Last peak length (nm) ; Max force (pN) ; Min force (pN) ; Max delta (nm) ; Min delta (nm) ; Peaks Diff\n')
90 f=open(self.my_work_dir+peackforce_filename,'w+')
91 f.write('Analysis started '+time.asctime()+'\n')
92 f.write('----------------------------------------\n')
93 f.write('; Peak number ; 1 peak Length (nm) ; 1 peak Force (pN) ; 2 peak Length (nm) ; 2 peak Force (pN) ; 3 peak Length (nm) ; 3 peak Force (pN) ; 4 peak Length (nm) ; 4 peak Force (pN) ; 5 peak Length (nm) ; 5 peak Force (pN) ; 6 peak Length (nm) ; 6 peak Force (pN) ; 7 peak Length (nm) ; 7 peak Force (pN) ; 8 peak Length (nm) ; 8 peak Force (pN)\n')
97 def plot_informations(itplot,pl_value):
100 contact_point.absolute_coords (2.4584142802103689e-007, -6.9647135616234017e-009)
101 peak_point.absolute_coords (3.6047748250571423e-008, -7.7142802788854212e-009)
102 other_fit_point.absolute_coords (4.1666139243838867e-008, -7.3759393477579707e-009)
103 peak_location [510, 610, 703, 810, 915, 1103]
104 peak_size [-1.2729111505202212e-009, -9.1632775347399312e-010, -8.1707438353929907e-010, -8.0335812578148904e-010, -8.7483955226387558e-010, -3.6269619757067322e-009]
105 params [2.2433999931959462e-007, 3.3230248825175678e-010]
106 fit_errors [6.5817195369767644e-010, 2.4415923138871498e-011]
108 fit_points=int(self.config['auto_fit_points']) # number of points to fit before the peak maximum <50>
110 T=self.config['temperature'] #temperature of the system in kelvins. By default it is 293 K. <301.0>
111 cindex=self.find_contact_point(itplot[0]) #Automatically find contact point <158, libhooke.ClickedPoint>
112 contact_point=self._clickize(itplot[0].vectors[1][0], itplot[0].vectors[1][1], cindex)
114 base_index_0=peak_location[-1]+self.fit_interval_nm(peak_location[-1], itplot[0], self.config['auto_right_baseline'],False)
115 self.basepoints.append(self._clickize(itplot[0].vectors[1][0],itplot[0].vectors[1][1],base_index_0))
116 base_index_1=self.basepoints[0].index+self.fit_interval_nm(self.basepoints[0].index, itplot[0], self.config['auto_left_baseline'],False)
117 self.basepoints.append(self._clickize(itplot[0].vectors[1][0],itplot[0].vectors[1][1],base_index_1))
118 self.basecurrent=self.current.path
119 boundaries=[self.basepoints[0].index, self.basepoints[1].index]
121 to_average=itplot[0].vectors[1][1][boundaries[0]:boundaries[1]] #y points to average
122 avg=np.mean(to_average)
123 return fit_points, contact_point, pl_value, T, cindex, avg
125 def features_peaks(itplot, peak, fit_points, contact_point, pl_value, T, cindex, avg):
127 calculate informations for each peak and add they in
128 c_lengths, p_lengths, sigma_c_lengths, sigma_p_lengths, forces, slopes
138 slope_span=int(self.config['auto_slope_span'])
140 peak_point=self._clickize(itplot[0].vectors[1][0],itplot[0].vectors[1][1],peak)
141 other_fit_point=self._clickize(itplot[0].vectors[1][0],itplot[0].vectors[1][1],peak-fit_points)
143 points=[contact_point, peak_point, other_fit_point]
145 params, yfit, xfit, fit_errors = self.wlc_fit(points, itplot[0].vectors[1][0], itplot[0].vectors[1][1], pl_value, T, return_errors=True)
148 delta_to_measure=itplot[0].vectors[1][1][peak-delta_force:peak+delta_force]
149 y=min(delta_to_measure)
151 slope=self.linefit_between(peak-slope_span,peak)[0]
152 #check fitted data and, if right, add peak to the measurement
153 if len(params)==1: #if we did choose 1-value fit
155 c_leng=params[0]*(1.0e+9)
157 sigma_c_leng=fit_errors[0]*(1.0e+9)
158 force = abs(y-avg)*(1.0e+12)
160 p_leng=params[1]*(1.0e+9)
161 #check if persistent length makes sense. otherwise, discard peak.
162 if p_leng>self.config['auto_min_p'] and p_leng<self.config['auto_max_p']:
164 p_lengths.append(p_leng)
165 c_lengths.append(params[0]*(1.0e+9))
166 sigma_c_lengths.append(fit_errors[0]*(1.0e+9))
167 sigma_p_lengths.append(fit_errors[1]*(1.0e+9))
168 forces.append(abs(y-avg)*(1.0e+12))
171 c_leng=params[0]*(1.0e+9)
172 sigma_c_leng=fit_errors[0]*(1.0e+9)
173 sigma_p_leng=fit_errors[1]*(1.0e+9)
174 force=abs(y-avg)*(1.0e+12)
178 #return c_lengths, p_lengths, sigma_c_lengths, sigma_p_lengths, forces, slopes
179 return c_leng, p_leng, sigma_c_leng, sigma_p_leng, force, slope
182 # ------ PROGRAM -------
184 for item in self.current_list:
186 item.identify(self.drivers)
187 itplot=item.curve.default_plots()
188 flatten=self._find_plotmanip('flatten') #extract flatten plot manipulator
189 itplot[0]=flatten(itplot[0], item, customvalue=1)
191 peak_location,peak_size=self.exec_has_peaks(item,min_deviation)
193 #We have troubles with exec_has_peaks (bad curve, whatever).
194 #Print info and go to next cycle.
195 print 'Cannot process ',item.path
198 if len(peak_location)==0:
202 fit_points, contact_point, pl_value, T, cindex, avg = plot_informations(itplot,pl_value)
204 print '\n\nCurve',item.path, 'is',c,'of',len(self.current_list),': found '+str(len(peak_location))+' peaks.'
206 #initialize output data vectors
214 #loop each peak of my curve
215 for peak in peak_location:
216 c_leng, p_leng, sigma_c_leng, sigma_p_leng, force, slope = features_peaks(itplot, peak, fit_points, contact_point, pl_value, T, cindex, avg)
217 for var, vector in zip([c_leng, p_leng, sigma_c_leng, sigma_p_leng, force, slope],[c_lengths, p_lengths, sigma_c_lengths, sigma_p_lengths, forces, slopes]):
221 #FIXME: We need a dictionary here...
222 allvects=[c_lengths, p_lengths, sigma_c_lengths, sigma_p_lengths, forces, slopes]
223 for vect in allvects:
225 for i in range(len(c_lengths)):
228 print 'Measurements for all peaks detected:'
229 print 'contour (nm)', c_lengths
230 print 'sigma contour (nm)',sigma_c_lengths
231 print 'p (nm)',p_lengths
232 print 'sigma p (nm)',sigma_p_lengths
233 print 'forces (pN)',forces
234 print 'slopes (N/m)',slopes
237 write automeasure text file
239 print 'Saving automatic measurement...'
240 f=open(self.my_work_dir+pclust_filename,'a+')
241 f.write(item.path+'\n')
242 for i in range(len(c_lengths)):
243 f.write(' ; '+str(c_lengths[i])+' ; '+str(p_lengths[i])+' ; '+str(forces[i])+' ; '+str(slopes[i])+' ; '+str(sigma_c_lengths[i])+' ; '+str(sigma_p_lengths[i])+'\n')
246 peak_number=len(c_lengths)
249 write peackforce text file
251 print 'Saving automatic measurement...'
252 f=open(self.my_work_dir+peackforce_filename,'a+')
253 f.write(item.path+'\n')
255 for i in range(len(c_lengths)):
256 peackforce_info = peackforce_info + ' ; ' + str(c_lengths[i]) + ' ; ' + str(forces[i])
257 f.write(' ; '+str(peak_number)+peackforce_info+'\n')
261 calculate clustering coordinates
265 for i in range(len(c_lengths)-1):
266 deltas.append(c_lengths[i+1]-c_lengths[i])
268 delta_mean=np.mean(deltas)
269 delta_median=np.median(deltas)
271 force_mean=np.mean(forces)
272 force_median=np.median(forces)
274 first_peak_cl=c_lengths[0]
275 last_peak_cl=c_lengths[-1]
277 max_force=max(forces[:-1])
278 min_force=min(forces)
280 max_delta=max(deltas)
281 min_delta=min(deltas)
283 delta_stdev=np.std(deltas)
284 forces_stdev=np.std(forces[:-1])
286 peaks_diff=(last_peak_cl-first_peak_cl)/peak_number
289 print 'Peaks',peak_number
290 print 'Mean delta',delta_mean
291 print 'Median delta',delta_median
292 print 'Mean force',force_mean
293 print 'Median force',force_median
294 print 'First peak',first_peak_cl
295 print 'Last peak',last_peak_cl
296 print 'Max force',max_force
297 print 'Min force',min_force
298 print 'Max delta',max_delta
299 print 'Min delta',min_delta
300 print 'Delta stdev',delta_stdev
301 print 'Forces stdev',forces_stdev
302 print 'Peaks difference',peaks_diff
305 write clustering coordinates
307 f=open(self.my_work_dir+realclust_filename,'a+')
308 f.write(item.path+'\n')
309 f.write(' ; '+str(peak_number)+ # non considerato
310 ' ; '+str(delta_mean)+ # 0
311 ' ; '+str(delta_median)+ # 1 -
312 ' ; '+str(force_mean)+ # 2
313 ' ; '+str(force_median)+ # 3 -
314 ' ; '+str(first_peak_cl)+ # 4 -
315 ' ; '+str(last_peak_cl)+ # 5 -
316 ' ; '+str(max_force)+ # 6
317 ' ; '+str(min_force)+ # 7
318 ' ; '+str(max_delta)+ # 8
319 ' ; '+str(min_delta)+ # 9
320 ' ; '+str(delta_stdev)+ # 10
321 ' ; '+str(forces_stdev)+ # 11
322 ' ; '+str(peaks_diff)+ # 12
327 self.do_pca(pclus_dir+"/"+realclust_filename)
330 def do_pca(self,args):
332 PCA -> "pca gaeta_coor_blind50.txt 1,3,6"
333 Automatically measures pca from coordinates filename and shows two interactives plots
334 With the second argument (arbitrary) you can select the columns and the multiplier factor
335 to use for the pca (for es "1,3*50,6,8x10,9"). Dont use spaces. "*" or "x" are the same thing.
336 Without second argument it reads pca_config.txt file
337 (c)Paolo Pancaldi, Massimo Sandal 2009
340 # reads the columns of pca
341 conf=open(config_file_path("pca_config.txt"), 'r')
342 config = conf.readlines()
345 self.plot_myCoord = [] # tiene le coordinate prese direttamente dal file creato con pCluster
346 self.plot_origCoord = [] # tiene le coordinate solo delle colonne scelte e moltiplicate per i valori scelti
347 self.plot_pcaCoord = [] # tiene le due colonne della pca
348 self.plot_pcaCoordTr = [] # tiene le due colonne della pca trasposta
349 self.plot_FiltOrigCoord = [] # tiene le coordinate solo dei punti filtrati per densita
350 self.plot_FiltPaths = [] # tiene i paths dei plot solo dei punti filtrati per densita
351 self.plot_paths = [] # tiene i paths di tutti i plots
352 self.plot_NewPcaCoord = [] # tiene le due colonne della pca filtrate per densita
353 self.plot_NewPcaCoordTr=[] # tiene le due colonne della pca trasposta filtrate per densita
356 # prende in inpunt un arg (nome del file)
357 # e il secondo le colonne su cui lavorare (e' arbitrario, riceve x es "1,2,3")
358 arg = args.split(" ")
365 # creo l'array "plot_myCoord" con tutte le coordinate dei plots
366 # e l'array plot_paths con tutti i percorsi dei plots
367 nPlotTot = -3 #tolgo le prime 3 righe iniziali del file
371 if row[0]!=" " and row[0]!="":
372 nPlotTot = nPlotTot+1
374 if row[0]==" " and row.find('nan')==-1 and row.find("-1.#IND")==-1:
375 row = row[row.index(";",2)+2:].split(" ; ") # non considero la prima colonna col #picchi
376 row = [float(i) for i in row]
378 #0:Mean delta, 1:Median delta, 2:Mean force, 3:Median force, 4:First peak length, 5:Last peak length
379 #6:Max delta 7:Min delta 8:Max force 9:Min force 10:Std delta 11:Std force
380 if (row[0]<500 and row[1]<500 and row[2]<500 and row[3]<500 and row[4]<500 and row[5]<500 and row[6]<500 and row[7]<500 and row[8]<500 and row[9]<500 and row[10]<500 and row[11]<500):
381 if (row[0]>0 and row[1]>0 and row[2]>0 and row[3]>0 and row[4]>0 and row[5]>0 and row[6]>0 and row[7]>0 and row[8]>0 and row[9]>0 and row[10]>0 and row[11]>0):
382 #row = row[0], row[2], row[3]*3, row[6], row[7]*56, row[8]
383 self.plot_myCoord.append(row)
384 self.plot_paths.append(plot_path_temp)
387 # creo l'array con alcune colonne e pure moltiplicate
388 for row in self.plot_myCoord:
390 for cols in config[0].split(","):
391 if cols.find("*")!=-1:
392 col = int(cols.split("*")[0])
393 molt = int(cols.split("*")[1])
394 elif cols.find("x")!=-1:
395 col = int(cols.split("x")[0])
396 molt = int(cols.split("x")[1])
400 res.append(row[col]*molt)
401 self.plot_origCoord.append(res)
403 # array convert, calculate PCA, transpose
404 self.plot_origCoord = np.array(self.plot_origCoord,dtype='float')
405 #print self.plot_origCoord.shape
406 self.plot_pcaCoord = pca(self.plot_origCoord, output_dim=2) #other way -> y = mdp.nodes.PCANode(output_dim=2)(array)
407 self.plot_pcaCoordTr = np.transpose(self.plot_pcaCoord)
408 pca_X=np.array(self.plot_pcaCoordTr[0],dtype='float')
409 pca_Y=np.array(self.plot_pcaCoordTr[1],dtype='float')
412 # Start section of testing with good plots # 4 TESTING!
419 goodnamefile=open(file_name.replace("coordinate", "good"),'r')
420 goodnames=goodnamefile.readlines()
421 nPlotGood = len(goodnames)-2 #tolgo prima e ultima riga
422 goodnames=[i.split()[0] for i in goodnames[1:]]
424 for index in range(len(self.plot_paths)):
425 if self.plot_paths[index][:-1] in goodnames:
426 Xsyn_1.append(pca_X[index])
427 Ysyn_1.append(pca_Y[index])
429 Xbad_1.append(pca_X[index])
430 Ybad_1.append(pca_Y[index])
431 # Stop section of testing with good plots # 4 TESTING!
435 clustplot1=lhc.PlotObject()
436 clustplot1.add_set(pca_X,pca_Y)
437 #clustplot1.add_set(Xbad_1,Ybad_1) # 4 TESTING!
438 #clustplot1.add_set(Xsyn_1,Ysyn_1) # 4 TESTING!
439 clustplot1.normalize_vectors()
440 clustplot1.styles=['scatter', 'scatter','scatter']
441 clustplot1.colors=[None,'red','green']
442 clustplot1.destination=0
443 self._send_plot([clustplot1])
444 self.clustplot1=clustplot1
446 # density and filer estimation
447 kernel = sp.stats.kde.gaussian_kde(sp.c_[pca_X,pca_Y].T)
449 for i in range(len(pca_X)):
450 kern_value = kernel.evaluate([pca_X[i],pca_Y[i]])
451 if tallest < kern_value:
452 tallest = float(kern_value)
453 if float(config[1]) == 0:
454 my_filter = float(tallest / 3.242311147)
456 my_filter = float(config[1])
458 # section useful only for graphic printing
463 mX, mY = sp.mgrid[xmin:xmax:100j, ymin:ymax:100j]
464 Z = sp.rot90(sp.fliplr(sp.reshape(kernel(sp.c_[mX.ravel(), mY.ravel()].T).T, mX.T.shape)))
465 axis_X = np.linspace(xmin,xmax,num=100)
466 axis_Y = np.linspace(ymin,ymax,num=100)
470 # tramite "kernel.evaluate" trovo lo score (altezza) di ogni coordinata e decido se mantenerla o no
473 filtered_PcaCoordTr = []
474 filtered_PcaCoord = []
475 for i in range(len(pca_X)):
476 kern_value = kernel.evaluate([pca_X[i],pca_Y[i]])
477 if kern_value > my_filter:
478 filtered_pca_X.append(pca_X[i])
479 filtered_pca_Y.append(pca_Y[i])
480 filtered_PcaCoordTr.append(filtered_pca_X)
481 filtered_PcaCoordTr.append(filtered_pca_Y)
482 filtered_PcaCoord = np.transpose(filtered_PcaCoordTr)
484 # creo i due array "plot_FiltOrigCoord" e "plot_FiltPaths" contenenti solo i dati filtrati con alta densita
485 for index in range(len(self.plot_pcaCoord)):
486 if self.plot_pcaCoord[index] in filtered_PcaCoord:
487 self.plot_FiltOrigCoord.append(self.plot_myCoord[index])
488 self.plot_FiltPaths.append(self.plot_paths[index])
491 # START PCA#2: USELESS!!!
493 # creo l array con alcune colonne e pure moltiplicate
495 for row in self.plot_FiltOrigCoord:
497 for cols in config[2].split(","):
498 if cols.find("*")!=-1:
499 col = int(cols.split("*")[0])
500 molt = int(cols.split("*")[1])
501 elif cols.find("x")!=-1:
502 col = int(cols.split("x")[0])
503 molt = int(cols.split("x")[1])
507 res.append(row[col]*molt)
508 temp_coord.append(res)
509 self.plot_FiltOrigCoord = temp_coord
511 # ricalcolo la PCA: array convert, calculate PCA, transpose
512 self.plot_FiltOrigCoord = np.array(self.plot_FiltOrigCoord,dtype='float')
513 #print self.plot_FiltOrigCoord.shape
514 self.plot_NewPcaCoord = pca(self.plot_FiltOrigCoord, output_dim=2) #other way -> y = mdp.nodes.PCANode(output_dim=2)(array)
515 self.plot_NewPcaCoordTr = np.transpose(self.plot_NewPcaCoord)
516 pca_X2=np.array(self.plot_NewPcaCoordTr[0],dtype='float')
517 pca_Y2=np.array(self.plot_NewPcaCoordTr[1],dtype='float')
519 # Start section of testing with good plots # 4 TESTING!
524 for index in range(len(self.plot_FiltPaths)):
525 if self.plot_FiltPaths[index][:-1] in goodnames:
526 Xsyn_2.append(pca_X2[index])
527 Ysyn_2.append(pca_Y2[index])
529 Xbad_2.append(pca_X2[index])
530 Ybad_2.append(pca_Y2[index])
533 clustplot2=lhc.PlotObject()
534 #clustplot2.add_set(pca_X2,pca_Y2)
535 clustplot2.add_set(Xbad_2,Ybad_2) # 4 TESTING!
536 clustplot2.add_set(Xsyn_2,Ysyn_2) # 4 TESTING!
537 clustplot2.normalize_vectors()
538 clustplot2.styles=['scatter', 'scatter','scatter']
539 clustplot2.colors=[None,'red','green']
540 clustplot2.destination=1
541 self._send_plot([clustplot2])
542 self.clustplot2=clustplot2
546 clustplot2=lhc.PlotObject()
547 clustplot2.add_set(filtered_pca_X,filtered_pca_Y)
548 clustplot2.normalize_vectors()
549 clustplot2.styles=['scatter', 'scatter','scatter']
550 clustplot2.colors=[None,'red','green']
551 clustplot2.destination=1
552 self._send_plot([clustplot2])
553 self.clustplot2=clustplot2
556 config_pca1 = config[0].replace("*", "x").rstrip("\n")
557 config_pca2 = config[2].replace("*", "x").rstrip("\n")
559 print "- START: "+file_name
560 print "Curve totali: ", nPlotTot
561 #print "Curve totali good: ", nPlotGood # 4 TESTING!
562 print "- FILTRO 1: 0-500 e NaN"
563 print "Curve totali rimaste: ", len(self.plot_origCoord)
564 #print 'Curve good rimaste: ', len(Xsyn_1) # 4 TESTING!
565 print "- FILTRO 2: PCA:"+config_pca1+" e DENSITA:"+str(my_filter)
566 print "Curve totali rimaste: ", len(self.plot_FiltOrigCoord)
567 #print 'Curve good rimaste: ', len(Xsyn_2) # 4 TESTING!
568 print "Piu alta: ", tallest
569 #print "- FILTRO 3: 2'PCA:"+config_pca2
572 # -- exporting coordinates and plot of PCA in debug mode! --
573 if config[3].find("true")!=-1:
574 #1' PCA: save plot and build coordinate s file
575 self.do_export(file_name.replace("coordinate_", "debug_pca1graph_").replace('.txt','_'+config_pca1) + " 0")
576 f = open(file_name.replace("coordinate_", "debug_pca1coor_").replace('.txt','_'+config_pca1+'.txt'),'w')
577 for i in range(len(pca_X)):
578 f.write (str(i) + "\t" + str(pca_X[i]) + "\t" + str(pca_Y[i]) + "\n")
580 #2' PCA: save plot and build coordinate s file
581 #self.do_export(file_name.replace("coordinate_", "debug_pca2graph_").replace('.txt','_'+config_pca2) + " 1")
582 #f = open(file_name.replace("coordinate_", "debug_pca2coor_").replace('.txt','_'+config_pca2+'.txt'),'w')
583 #for i in range(len(pca_X2)):
584 # f.write (str(i) + "\t" + str(pca_X2[i]) + "\t" + str(pca_Y2[i]) + "\n")
587 self.do_export(file_name.replace("coordinate_", "debug_densitygraph_").replace('.txt','_'+config_pca1+'_'+str(my_filter).replace(".",",")) + " 1")
588 f = open(file_name.replace("coordinate_", "debug_densitycoor_").replace('.txt','_'+config_pca1+'_'+str(my_filter).replace(".",",")+'.txt'),'w')
589 for i in range(len(filtered_pca_X)):
590 f.write (str(i) + "\t" + str(filtered_pca_X[i]) + "\t" + str(filtered_pca_Y[i]) + "\n")
592 #ALL GOOD COORDINATES (without NaN and 0<x<500)
593 f = open(file_name.replace("coordinate_", "debug_allgoodcoor_"),'w')
594 for i in range(len(self.plot_myCoord)):
595 for cel in self.plot_myCoord[i]:
596 f.write (" ; " + str(cel))
602 pcl_name = file_name.replace("coordinate_", "goodplots_").replace('.txt','_'+config_pca1+'_'+str(my_filter).replace(".",","))
603 if os.path.exists(pcl_name+slash): shutil.rmtree(pcl_name)
604 os.mkdir(pcl_name+slash)
605 f = open(pcl_name+'.txt','w')
606 for i in range(len(self.plot_FiltPaths)):
607 myfile = str(self.plot_FiltPaths[i]).rstrip("\n")
608 f.write (myfile+"\n")
609 shutil.copy2(myfile, pcl_name)
613 def do_multipca(self,args):
615 MULTIPCA -> "multipca gaeta_coor_blind50.txt 3"
616 Automatically multiply the column suggest in second argument for value between 1-100 (step of 2),
617 measures pca from coordinates filename and save the png plots.
618 (c)Paolo Pancaldi, Massimo Sandal 2009
620 # reads the columns of pca
621 conf=open("pca_config.txt")
622 config = conf.readlines() # config[0] = "1,2,3"
625 arg = args.split(" ")
628 for i in range(1, 51, 1):
629 self.do_pca(file_name + " " + config[0].replace(column,column+"*"+str(i),1))
631 def do_doublepca(self,args):
633 DOUBLEPCA -> "doublepca gaeta_coor_blind50.txt"
634 Automatically it launches the pca command for all combinations with two column
635 (c)Paolo Pancaldi, Massimo Sandal 2009
638 arg = args.split(" ")
640 for i in range(1, 13):
641 for j in range(1, 13):
643 self.do_pca(file_name + " " + str(i) + "," + str(j))
645 def do_triplepca(self,args):
647 TRIPLEPCA -> "triplepca gaeta_coor_blind50.txt"
648 Automatically it launches the pca command for all combinations with three column
649 (c)Paolo Pancaldi, Massimo Sandal 2009
652 arg = args.split(" ")
654 for i in range(1, 13):
655 for j in range(1, 13):
656 for k in range(1, 13):
657 if i!=j and i!=k and j!=k:
658 self.do_pca(file_name + " " + str(i) + "," + str(j) + "," + str(k))
660 def do_pclick(self,args):
662 It returns id, coordinates and file name of a clicked dot on a PCA graphic
665 self._send_plot([self.clustplot1]) #quick workaround for BAD problems in the GUI
667 point = self._measure_N_points(N=1, whatset=0)
668 indice = point[0].index
669 plot_file = self.plot_paths[indice]
670 dot_coord = self.plot_pcaCoord[indice]
671 print "file: " + str(plot_file).rstrip()
672 print "id: " + str(indice)
673 print "coord: " + str(dot_coord)
674 self.do_genlist(str(plot_file))
675 #self.do_jump(str(plot_file))
677 # indea iniziata e messa da parte...
678 def do_peakforce(self, args):
680 peackforce -> "peackforce peackforce_file.txt"
681 Automatically measures peack and force plots
682 (c)Paolo Pancaldi, Massimo Sandal 2009
685 # prende in inpunt un arg (nome del file)
689 # scrivo un file temp
690 g = open('_prove.txt','w')
697 if row[0]==" " and row.find('nan')==-1 and row.find("-1.#IND")==-1:
698 # FILTRO SUI 7 PICCHI
699 num_pic = int(row.split(" ; ")[1])
701 width_force = row.split(" ; ")
702 w1 = float(width_force[2]); f1 = float(width_force[3]);
703 w2 = float(width_force[4]); f2 = float(width_force[5]);
704 w3 = float(width_force[6]); f3 = float(width_force[7]);
705 w4 = float(width_force[8]); f4 = float(width_force[9]);
706 w5 = float(width_force[10]); f5 = float(width_force[11]);
707 w6 = float(width_force[12]); f6 = float(width_force[13]);
708 w7 = float(width_force[14]); f7 = float(width_force[15]);
709 if w1>0 and w1<1000 and w2>0 and w2<1000 and w3>0 and w3<1000 and w4>0 and w4<1000 and w5>0 and w5<1000 and w6>0 and w6<1000 and w7>0 and w7<1000:
710 score_76 = abs(32 - (w7 - w6))
711 score_65 = abs(32 - (w6 - w5))
712 score_54 = abs(32 - (w5 - w4))
713 score_43 = abs(32 - (w4 - w3))
714 score_32 = abs(32 - (w3 - w2))
715 score_21 = abs(32 - (w2 - w1))
716 writeme = str(score_76) + " --- " + str(row)