Various adjustments. I should commit more often ;).

Added # comments to tweakfile syntax.

Played around with adding a white-noise floor in the vibration fitting,
but the fits didn't look all that convincing.  Some of the white-noise
code is still in place, but I think it's currently disabled ;).

Fixed some typos in TEXT_VERBOSE output in vib_analyze.py

diff --git a/calibcant/T_analyze.py b/calibcant/T_analyze.py
index 6f0a8a12795dac041b2ba2108e1b9b0383adc479..d4ddf726f07b7aa07a0082eeba20f37d25860509 100755 (executable)
--- a/calibcant/T_analyze.py
+++ b/calibcant/T_analyze.py
@@ -106,7 +106,9 @@ def T_load_analyze_tweaked(tweak_file, convert_to_K=C_to_K, textVerboseFile=None
     Ts = []
     for line in file(tweak_file, 'r') :
         parsed = line.split()
-        path = parsed[0].split('\n')[0]
+        path = parsed[0].strip()
+        if path[0] == '#': # a comment
+            continue
         if textVerboseFile != None :
             print >> textVerboseFile, "Reading data from %s" % (path)
         # read the data
@@ -142,6 +144,7 @@ if __name__ == '__main__' :
                     'Tweak file mode:\n'
                     'Runs the same analysis as in single file mode for each T file in\n'
                     'a tweak file.  Each line in the tweak file specifies a single T file.\n'
+                    'Blank lines and those beginning with a pound sign (#) are ignored.\n'
                     'A T file contains a sequence of 32 bit floats representing temperature in K.\n'
                     )
     parser = OptionParser(usage=usage_string, version='%prog 0.1')

diff --git a/calibcant/analyze.py b/calibcant/analyze.py
index 7ed75645f0b3d0b1028a9b4dfbf24f574ba1938d..c4850095366319d370845b157cc1c35fe8e47e0a 100755 (executable)
--- a/calibcant/analyze.py
+++ b/calibcant/analyze.py
@@ -269,19 +269,24 @@ if __name__ == '__main__' :
 
     options,args = parser.parse_args()
     parser.destroy()
-
+    
+    config.TEXT_VERBOSE = options.verbose
+    config.PYLAB_INTERACTIVE = False
+    config.PYLAB_VERBOSE = options.plot
+    
     bumps = get_array(options.bump_string, options.bump_file, "bump")
     temps = get_array(options.temp_string, options.temp_file, "temp")
     vibs = get_array(options.vib_string, options.vib_file, "vib")
 
-    if options.plot :
-        calib_plot(bumps, temps, vibs)
+    #if options.plot :
+    #    calib_plot(bumps, temps, vibs)
 
     if options.celsius :
         for i in range(len(temps)) :
             temps[i] = T_analyze.C_to_K(temps[i])
 
-    km,ks,ps2m,ps2s,Tm,Ts,ooVp2m,ooVp2s = calib_analyze(bumps, temps, vibs)
+    km,ks,ps2m,ps2s,Tm,Ts,ooVp2m,ooVp2s = \
+        calib_analyze(bumps, temps, vibs, plotVerbose=options.plot)
 
     if options.verbose :
         print >> sys.stderr, \

diff --git a/calibcant/bump_analyze.py b/calibcant/bump_analyze.py
index 11a8f51e719f7f8b96220bec7fca47b4357ddb06..36a6b81d20edb5b287fdb1e25801ae960dc9b1fc 100755 (executable)
--- a/calibcant/bump_analyze.py
+++ b/calibcant/bump_analyze.py
@@ -128,7 +128,9 @@ def bump_load_analyze_tweaked(tweak_file, **kwargs):
     photoSensitivity = []
     for line in file(tweak_file, 'r') :
         parsed = line.split()
-        path = parsed[0].split('\n')[0]
+        path = parsed[0].strip()
+        if path[0] == '#' : # a comment
+            continue
         if config.TEXT_VERBOSE :
             print "Reading data from %s with ranges %s" % (path, parsed[1:])
         # read the data
@@ -192,6 +194,7 @@ if __name__ == '__main__' :
                     'Tweak file mode:\n'
                     'Runs the same analysis as in single file mode for each bump in\n'
                     'a tweak file.  Each line in the tweak file specifies a single bump.\n'
+                    'Blank lines and those beginning with a pound sign (#) are ignored.\n'
                     'The format of a line is a series of whitespace-separated fields--\n'
                     'a base file path followed by optional point index ranges, e.g.:\n'
                     '20080919/20080919132500_bump_surface 10:651 1413:2047\n'

diff --git a/calibcant/common.py b/calibcant/common.py
index c58e06b9a1208a627e796854d69bc9da92632c5e..1e8ff7137e3c839f6a4ed3496842428099bc59d5 100644 (file)
--- a/calibcant/common.py
+++ b/calibcant/common.py
@@ -31,7 +31,7 @@ _final_flush_plot = None
 _pylab = None
 def _dummy_fn(): pass
 
-def _import_pylab() :
+def _import_pylab() :  # TODO: auto detect no DISPLAY and abort
     """Import pylab plotting functions for when we need to plot.  This
     function can be called multiple times, and ensures that the pylab
     setup is only imported once.  It defines the functions

diff --git a/calibcant/config.py b/calibcant/config.py
index 8d7d51131a525536e4500d4003db17f62b1dd623..4dbd1a20b924ff4ca81a63bd2937d0441544cdb7 100644 (file)
--- a/calibcant/config.py
+++ b/calibcant/config.py
@@ -31,8 +31,8 @@ LOG_DATA = True  # quietly grab all real-world data and log to LOG_DIR
 LOG_DIR = '$DEFAULT$/calibrate_cantilever'
 GNUFIT_DATA_BASE='./calibrate_cantilever_fitdata'
 TEXT_VERBOSE = True      # for debugging
-GNUPLOT_VERBOSE = True     # turn on fit check plotting
-PYLAB_VERBOSE = True     # turn on plotting
+GNUPLOT_VERBOSE = True   # turn on fit check plotting
+PYLAB_VERBOSE = False    # turn on plotting
 PYLAB_INTERACTIVE = True # select between draw() and show() for flushing plots
 BASE_FIGNUM = 20 # to avoid writing to already existing figures
 

diff --git a/calibcant/vib_analyze.py b/calibcant/vib_analyze.py
index a98ea0534b9008e129d4176d1158487da8351b22..ff37f1242bed7b381463a3e690d822780f9ede1c 100755 (executable)
--- a/calibcant/vib_analyze.py
+++ b/calibcant/vib_analyze.py
@@ -103,13 +103,13 @@ def vib_analyze(deflection_bits, freq,
         FFT_tools.unitary_avg_power_spectrum((Vphoto_t - Vphoto_t.mean())*1e6,
                                              freq, chunk_size,overlap, window)
 
-    A,B,C = fit_psd(freq_axis, power, **kwargs)
+    A,B,C,D = fit_psd(freq_axis, power, **kwargs)
 
     if config.TEXT_VERBOSE : 
         print "Fitted f(x) = C / ((A**2-x**2)**2 + (x*B)**2) with"
-        print "A = %g, \t B = %g, \t C = %g" % (A, B, C)
+        print "A = %g, \t B = %g, \t C = %g, \t D = %g" % (A, B, C, D)
 
-    vib_plot(deflection_bits, freq_axis, power, A, B, C, **kwargs)
+    vib_plot(deflection_bits, freq_axis, power, A, B, C, D, **kwargs)
 
     # Our A is in uV**2, so convert back to Volts**2
     fitted_variance = lorentzian_area(A,B,C) * 1e-12
@@ -117,8 +117,8 @@ def vib_analyze(deflection_bits, freq,
     naive_variance = vib_analyze_naive(deflection_bits,
                                        Vphoto_in2V=Vphoto_in2V)
     if config.TEXT_VERBOSE:
-        print "Fitted variance: %g V**2", fitted_variance
-        print "Naive variance : %g V**2", naive_variance
+        print "Fitted variance: %g V**2" % fitted_variance
+        print "Naive variance : %g V**2" % naive_variance
     
     return min(fitted_variance, naive_variance)
 
@@ -175,6 +175,7 @@ def fit_psd(freq_axis, psd_data, minFreq=500, maxFreq=7000) :
     B_guess = res_freq / numpy.sqrt(Q_guess**2-0.5)
     A_guess = Q_guess*B_guess
     C_guess = max_psd * (-res_freq**4 + A_guess**4)
+    D_guess = 0 # allow a noise floor offset (DISABLED in fitting below)
     # 
     # Half width w on lower side when L(a-w) = L(a)/2
     #  (a**2 - (a-w)**2)**2 + (b*(a-w))**2 = 2*(b*a)**2
@@ -241,9 +242,10 @@ def fit_psd(freq_axis, psd_data, minFreq=500, maxFreq=7000) :
             else :
                 print "Solution did not converge"
         A,B,C = p
+        D = 0 # do not fit the noise floor (fit doesn't look very convincing)
         A=abs(A) # A and B only show up as squares in f(x)
         B=abs(B) # so ensure we get positive values
-    return (A, B, C)
+    return (A, B, C, D)
 
 def lorentzian_area(A, B, C) :
     # Integrating the the power spectral density per unit time (power) over the
@@ -290,7 +292,7 @@ def vib_load(datafile=None) :
     return data
 
 @splittableKwargsFunction()
-def vib_plot(deflection_bits, freq_axis, power, A, B, C,
+def vib_plot(deflection_bits, freq_axis, power, A, B, C, D,
              minFreq=None, maxFreq=None, plotVerbose=False) :
     """
     If plotVerbose or config.PYLAB_VERBOSE == True, plot 3 subfigures:
@@ -342,8 +344,10 @@ def vib_plot(deflection_bits, freq_axis, power, A, B, C,
             common._pylab.axvspan(minFreq, maxFreq, facecolor='g', alpha=0.1,
                                   zorder = -1)
         
-        fitdata = C / ((A**2-freq_axis**2)**2 + (B*freq_axis)**2)
-        common._pylab.plot(freq_axis, fitdata, 'b-');
+        fitdata = C / ((A**2-freq_axis**2)**2 + (B*freq_axis)**2) + D
+        common._pylab.plot(freq_axis, fitdata, 'b-')
+        noisefloor = D + 0*freq_axis;
+        common._pylab.plot(freq_axis, noisefloor)
         common._pylab.hold(False)
         axes.axis([xmin,xmax,ymin,ymax])
 
@@ -381,7 +385,7 @@ def vib_plot(deflection_bits, freq_axis, power, A, B, C,
 # were defined).
 make_splittable_kwargs_function(vib_analyze,
     (fit_psd, 'freq_axis', 'psd_data'),
-    (vib_plot, 'deflection_bits', 'freq_axis', 'power', 'A', 'B', 'C',
+    (vib_plot, 'deflection_bits', 'freq_axis', 'power', 'A', 'B', 'C', 'D',
      'minFreq', 'maxFreq'))
 
 @splittableKwargsFunction((vib_analyze, 'deflection_bits', 'freq'))
@@ -395,8 +399,9 @@ def vib_load_analyze_tweaked(tweak_file, **kwargs) :
     dl = data_logger.data_load()
     Vphoto_var = []
     for path in file(tweak_file, 'r') :
-        if path[-1] == '\n':
-            path = path.split('\n')[0]
+        path = path.strip()
+        if path[0] == '#': # a comment
+            continue
         # read the data
         data = vib_load(path)
         deflection_bits = data['Deflection input']
@@ -432,7 +437,9 @@ if __name__ == '__main__' :
                     'Deflection power spectral density (PSD) data (X^2/Hz)\n'
                     'and returns the variance in X units (X^2)\n'
                     '<input-file> should be 1 column ASCII without a header line.\n'
-                    'e.g: "<deflection bits>\\n..."\n')
+                    'e.g: "<deflection bits>\\n..."\n'
+                    #TODO, better overview of input modes, e.g. a la T_analyze
+                    )
     parser = OptionParser(usage=usage_string, version='%prog 0.1')
     parser.add_option('-f', '--sample-freq', dest='freq', default=100e3,
                       help='sample frequency in Hz (default %default)',