from numpy import logspace, log10, mean, linspace, ceil
from connectClient import CEE
import pylab
import time

CEE = CEE()
pylab.ion()

periodCount = 10 

sampleTime = CEE.devInfo['sampleTime']
minFreq = 10
maxFreq = (1/sampleTime)/20
# maximum frequency = sampling rate / 20

frequencies = logspace(log10(minFreq), log10(maxFreq), 40)
#log-spaced array of frequencies
amplitudes = []
phases = []

def findPhases(maxesA, maxesB):
	phases = []
	maxesA = maxesA[1::]
	maxesB = maxesB[1::]
	for i in range(len(maxesA)):
		phase = ( maxesA[i] - maxesB[i] ) * frequency * 360.0
		if phase > 180:
			phase = 180 - phase
		if phase < -180:
			phase = 360 + phase
		phases.append(phase)
	return phases 

def findLocalMaxes(values, a):
	chunkwidth = len(values)/periodCount
	# determine samples per period
	split = [values[i:i+chunkwidth] for i in range(0, len(values), chunkwidth)]
	# split into periods
	localMaxes = map(max, split)
	# find local maximums
	localMaxTimes = [ (chunk.index(localMax) + split.index(chunk)*chunkwidth) * sampleTime
		for localMax, chunk in zip(localMaxes, split) ]
	# find indexes of local maximums
	return localMaxes, localMaxTimes
	# trim the garbage datapoints and return the information

for frequency in frequencies:
	period = 1/frequency
	setResponse = CEE.setOutput('a', 'v', 2.5, 'sine', 2.5, frequency, 0, 0)
	# source sine wave with full-scale voltage range at target frequency
	sampleCount = ( period * periodCount ) / sampleTime
	# do math to get the equivalent of 'periodCount' in samples
	startSample = ceil(setResponse['startSample']/period)*period
	v = CEE.getInput('a', 0, int(sampleCount), int(startSample))[0]
	i = CEE.getInput('b', 0, int(sampleCount), int(startSample))[0]
	# get samples from CEE
	vMaxes, vMaxTimes = findLocalMaxes(v, True)
	iMaxes, iMaxTimes = findLocalMaxes(i, False)
	# get timestamps and maximum values for voltage and current
	amplitudes.append(mean(iMaxes))
	# calculate and record amplitude from array of local maximum currents
	phases.append(mean(findPhases(vMaxTimes, iMaxTimes)))
	# calculate and record phases from v/i maximums' timestamps 

pylab.figure()
pylab.subplot(2,1,1)
pylab.semilogx(frequencies, amplitudes, '.')
pylab.xlim(minFreq,maxFreq)
pylab.ylim(0,5)
pylab.ylabel("mean peak voltage")
pylab.subplot(2,1,2)
pylab.semilogx(frequencies, phases, '.')
pylab.xlim(minFreq,maxFreq)
pylab.ylabel("phase shift in degrees")
pylab.xlabel("frequency")

@@ -36,24 +36,8 @@
 mean
-, linspace, ceil
 %0Afro
@@ -314,13 +314,15 @@
 q), 
-4
+10
 0)%0A#
+ 
 log-
@@ -427,52 +427,8 @@
  %5B%5D%0A
-%09maxesA = maxesA%5B1::%5D%0A%09maxesB = maxesB%5B1::%5D%0A
 %09for
@@ -490,16 +490,29 @@
 B%5Bi%5D ) *
+ sampleTime *
  frequen
@@ -522,16 +522,36 @@
 * 360.0%0A
+%09%09# calculate phase%0A
 %09%09if pha
@@ -625,16 +625,36 @@
 + phase%0A
+%09%09# normalize phase%0A
 %09%09phases
@@ -710,19 +710,16 @@
 s(values
-, a
 ):%0A%09chun
@@ -1009,29 +1009,16 @@
 nkwidth)
- * sampleTime
 %0A%09%09for l
@@ -1112,16 +1112,21 @@
 calMaxes
+%5B1::%5D
 , localM
@@ -1132,16 +1132,21 @@
 MaxTimes
+%5B1::%5D
 %0A%09# trim
@@ -1511,70 +1511,8 @@
 les%0A
-%09startSample = ceil(setResponse%5B'startSample'%5D/period)*period%0A
 %09v =
@@ -1543,36 +1543,45 @@
 t(sampleCount), 
-int(
+setResponse%5B'
 startSample))%5B0%5D
@@ -1571,25 +1571,26 @@
 'startSample
-)
+'%5D
 )%5B0%5D%0A%09i = CE
@@ -1626,20 +1626,29 @@
 Count), 
-int(
+setResponse%5B'
 startSam
@@ -1650,17 +1650,18 @@
 rtSample
-)
+'%5D
 )%5B0%5D%0A%09# 
@@ -1722,14 +1722,8 @@
 es(v
-, True
 )%0A%09i
@@ -1761,15 +1761,8 @@
 es(i
-, False
 )%0A%09#

clean and comment