hooke/plugin/jumpstat.py

   1 # Copyright
   2
   3 from ..libhooke import WX_GOOD, ClickedPoint
   4 import wxversion
   5 wxversion.select(WX_GOOD)
   6 from wx import PostEvent
   7 import numpy as np
   8 import scipy as sp
   9 import copy
  10 import os.path
  11 import time
  12 import sys
  13 import warnings
  14 warnings.simplefilter('ignore',np.RankWarning)
  15
  16
  17 class jumpstatCommands():
  18
  19     def do_jumpstat(self,args):
  20      '''
  21      JUMPSTAT
  22      jumpstat.py
  23      Based on the convolution recognition automatically give:
  24      - the delta distance between the peaks,
  25      - the delta-force from the top of the peaks and subsequent relaxation,
  26      - the delta-force from the top of the peaks and the baseline
  27      The command allow also to remove the unwanted peaks that can be due to interference.
  28      When you first issue the command, it will ask for the filename. If you are giving the filename
  29      of an existing file, autopeak will resume it and append measurements to it. If you are giving
  30      a new filename, it will create the file and append to it until you close Hooke.
  31      You can also define a minimun deviation of the peaks.
  32
  33      Syntax:
  34      jumpstat [deviation]
  35      deviation = number of times the convolution signal is above the noise absolute deviation.
  36      '''
  37
  38
  39      #finding the max and the minimum positions for all the peaks
  40      noflatten=False
  41      #we use if else only to avoid a "bad input" message from find_current_peaks
  42      if (len(args)==0):
  43              max_peaks_location, peak_size=self.find_current_peaks(noflatten)
  44              min_peaks_location, pks2=self.find_current_peaks(noflatten, True, False)
  45      else:
  46              max_peaks_location, peak_size=self.find_current_peaks(noflatten, args)
  47              min_peaks_location, pks2=self.find_current_peaks(noflatten, args, False)
  48
  49
  50      #print "max_peaks_location:  "+str(len(max_peaks_location))
  51      #print "min_peaks_location:  "+str(len(min_peaks_location))
  52
  53      #if no peaks, we have nothing to plot. exit.
  54      if len(max_peaks_location)==0:
  55             print "No peaks on this curve."
  56             return
  57
  58      if len(max_peaks_location)!=len(min_peaks_location):
  59             print "Something went wrong in peaks recognition, number of minima is different from number of maxima. Exiting."
  60             return
  61
  62      #otherwise, we plot the peak locations.
  63      xplotted_ret=self.plots[0].vectors[1][0]
  64      yplotted_ret=self.plots[0].vectors[1][1]
  65      xgood=[xplotted_ret[index] for index in max_peaks_location]
  66      ygood=[yplotted_ret[index] for index in max_peaks_location]
  67
  68      xafter=[xplotted_ret[index] for index in min_peaks_location]
  69      yafter=[yplotted_ret[index] for index in min_peaks_location]
  70
  71      recplot=self._get_displayed_plot()
  72      recplot2=self._get_displayed_plot()
  73      recplot.vectors.append([xgood,ygood])
  74      recplot2.vectors.append([xafter,yafter])
  75
  76      if recplot.styles==[]:
  77          recplot.styles=[None,None,'scatter']
  78          recplot.colors=[None,None,None]
  79      else:
  80          recplot.styles+=['scatter']
  81          recplot.colors+=[None]
  82
  83      if recplot2.styles==[]:
  84          recplot2.styles=[None,None,None]
  85          recplot2.colors=[None,'1.0',None]
  86      else:
  87          recplot2.styles+=['scatter']
  88          recplot2.colors+=['0.5']
  89
  90      self._send_plot([recplot])
  91      self._send_plot([recplot2])
  92
  93
  94      #finding the baseline
  95      self.basepoints=self.baseline_points(max_peaks_location, recplot)
  96      boundaries=[self.basepoints[0].index, self.basepoints[1].index]
  97      boundaries.sort()
  98      to_average=recplot.vectors[1][1][boundaries[0]:boundaries[1]] #y points to average
  99      avg=np.mean(to_average)
 100
 101
 102      dist=[]
 103      jumpforce=[]
 104      force=[]
 105
 106      #we calculate the distance vector
 107      for g in range(len(max_peaks_location)-1):
 108        dist.append((10**9)*(xplotted_ret[max_peaks_location[g]]-xplotted_ret[max_peaks_location[g+1]]))
 109      print "Distance values for the peaks in nm:"
 110      print dist
 111
 112      #the jump-force vector
 113      for g in range(len(max_peaks_location)):
 114       jumpforce.append((10**12)     *(yplotted_ret[min_peaks_location[g]] -yplotted_ret[max_peaks_location[g]])   )
 115      print "Force values for the jumps of the peaks in pN:"
 116      print jumpforce
 117
 118      #the force from baseline vector
 119      for g in range(len(max_peaks_location)):
 120       force.append((10**12)*(avg-yplotted_ret[max_peaks_location[g]]))
 121      print "Force values for the peaks in pN:"
 122      print force
 123
 124
 125
 126      #Now ask for the peaks that we don't want
 127      print 'Peaks to ignore (0,1...n from contact point,return to take all)'
 128      print 'N to discard measurement'
 129      exclude_raw=raw_input('Input:')
 130      if exclude_raw=='N':
 131         print 'Discarded.'
 132         return
 133
 134      if not exclude_raw=='':
 135         exclude=exclude_raw.split(',')
 136         #we convert in numbers the input
 137         try:
 138             exclude=[int(item) for item in exclude]
 139         except:
 140             print 'Bad input, taking nothing.'
 141             return
 142
 143 #    we remove the peaks that we don't want from the list, we need a counter beacause if we remove
 144 #    a peaks the other peaks in the list are shifted by one at each step
 145         count=0
 146         for a in exclude:
 147           if (a==0):
 148              max_peaks_location=max_peaks_location[1:]
 149              min_peaks_location=min_peaks_location[1:]
 150           else:
 151              new_a=a-count
 152              max_peaks_location=  max_peaks_location[0:new_a]+max_peaks_location[new_a+1:]
 153              min_peaks_location=  min_peaks_location[0:new_a]+min_peaks_location[new_a+1:]
 154              peak_size=            peak_size[0:new_a]+peak_size[new_a+1:]
 155           count+=1
 156
 157
 158      #print "max_peaks_location:  "+str(len(max_peaks_location))
 159      #print "min_peaks_location:  "+str(len(min_peaks_location))
 160
 161
 162      dist=[]
 163      jumpforce=[]
 164      force=[]
 165      #we recalculate the distances and the forces after the removing of the unwanted peaks
 166      for g in range(len(max_peaks_location)-1):
 167          dist.append(xplotted_ret[max_peaks_location[g]]-xplotted_ret[max_peaks_location[g+1]])
 168      for g in range(len(max_peaks_location)):
 169          jumpforce.append( yplotted_ret[min_peaks_location[g]] - yplotted_ret[max_peaks_location[g]]   )
 170      for g in range(len(max_peaks_location)):
 171          force.append(avg  -  yplotted_ret[max_peaks_location[g]])
 172
 173
 174
 175
 176
 177         #Save file info
 178      if self.autofile=='':
 179             self.autofile=raw_input('Jumpstat filename? (return to ignore) ')
 180             if self.autofile=='':
 181                 print 'Not saved.'
 182                 return
 183
 184      if not os.path.exists(self.autofile):
 185             f=open(self.autofile,'w+')
 186             f.write('Analysis started '+time.asctime()+'\n')
 187             f.write('----------------------------------------\n')
 188             f.write('; Delta Distance length (m);  Jump Force pN;  Standard Force pN\n')
 189             f.write(self.current.path+'\n')
 190             for k in range(len(dist)):
 191                f.write(";")
 192                f.write(str(dist[k])+";"+str(jumpforce[k])+";"+str(force[k])+"\n"   )
 193             f.write("\n")
 194             f.close()
 195
 196      else:
 197             f=open(self.autofile,'a+')
 198             f.write(self.current.path+'\n')
 199             for k in range(len(dist)):
 200               f.write(";")
 201               f.write(str(dist[k])+";"+str(jumpforce[k])+";"+str(force[k])+"\n"   )
 202             f.write("\n")
 203             f.close()
 204
 205      print 'Saving...'