hooke/plugin/jumpstat.py

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