7 if len(sys.argv) != 2 and len(sys.argv) != 4:
8 print("Usage : %s datafile", sys.argv[0])
9 print("or : %s datafile p1 p2", sys.argv[0])
10 print("where : p1 < p2 belongs to sizes in datafiles")
13 if len(sys.argv) == 4:
17 ##-----------------------------------------
18 ## avg : return average of a list of values
19 ## param l list of values
20 ##-----------------------------------------
27 ##-------------------------------------------------
29 ## param X first data vector (..x_i..)
30 ## param Y second data vector (..x_i..)
31 ## = 1/n \Sum_{i=1}^n (x_i - avg(x)) * (y_i - avg(y))
32 ##--------------------------------------------------
35 n=len(X) # n=len(X)=len(Y)
40 S_XY = S_XY + ((X[i]-avg_X)*(Y[i]-avg_Y))
45 ##----------------------------------
46 ## variance : variance
47 ## param X data vector ( ..x_i.. )
48 ## (S_X)^2 = (Sum ( x_i - avg(x) )^2 ) / n
49 ##----------------------------------
55 S_X2 = S_X2 + ((X[i] - avg_X)**2)
59 ##-----------------------------------------------------------------------------------------------
60 ## correl_split_weighted : compute regression on each segment and
61 ## return the weigthed sum of correlation coefficients
62 ## param X first data vector (..x_i..)
63 ## param Y second data vector (..x_i..)
64 ## param segments list of pairs (i,j) where i refers to the ith value in X, and jth value in X
65 ## return (C,[(i1,j1,X[i1],X[j1]), (i2,j2,X[i2],X[j2]), ....]
66 ## where i1,j1 is the first segment, c1 the correlation coef on this segment, n1 the number of values
67 ## i2,j2 is the second segment, c2 the correlation coef on this segment, n2 the number of values
69 ## and C=c1/n1+c2/n2+...
70 ##-----------------------------------------------------------------------------------------------
71 def correl_split_weighted( X , Y , segments ):
72 # expects segments = [(0,i1-1),(i1-1,i2-1),(i2,len-1)]
74 interv = list(); # regr. line coeffs and range
77 for (start,stop) in segments:
78 sum_nb_val = sum_nb_val + stop - start;
81 S_XY= cov( X [start:stop+1], Y [start:stop+1] )
82 S_X2 = variance( X [start:stop+1] )
83 S_Y2 = variance( Y [start:stop+1] ) # to compute correlation
86 c = S_XY/(sqrt(S_X2)*sqrt(S_Y2))
87 a = S_XY/S_X2 # regr line coeffs
88 b= avg ( Y[start:stop+1] ) - a * avg( X[start:stop+1] )
89 print(" range [%d,%d] corr=%f, coeff det=%f [a=%f, b=%f]" % (X[start],X[stop],c,c**2,a, b))
90 correl.append( (c, stop-start) ); # store correl. coef + number of values (segment length)
91 interv.append( (a,b, X[start],X[stop]) );
94 glob_corr = glob_corr + (l/sum_nb_val)*c # weighted product of correlation
95 print('-- %f * %f' % (c,l/sum_nb_val))
97 print("-> glob_corr=%f\n" % glob_corr)
98 return (glob_corr,interv);
103 ##-----------------------------------------------------------------------------------------------
104 ## correl_split : compute regression on each segment and
105 ## return the product of correlation coefficient
106 ## param X first data vector (..x_i..)
107 ## param Y second data vector (..x_i..)
108 ## param segments list of pairs (i,j) where i refers to the ith value in X, and jth value in X
109 ## return (C,[(i1,j1,X[i1],X[j1]), (i2,j2,X[i2],X[j2]), ....]
110 ## where i1,j1 is the first segment, c1 the correlation coef on this segment,
111 ## i2,j2 is the second segment, c2 the correlation coef on this segment,
114 ##-----------------------------------------------------------------------------------------------
115 def correl_split( X , Y , segments ):
116 # expects segments = [(0,i1-1),(i1-1,i2-1),(i2,len-1)]
118 interv = list(); # regr. line coeffs and range
120 for (start,stop) in segments:
123 S_XY= cov( X [start:stop+1], Y [start:stop+1] )
124 S_X2 = variance( X [start:stop+1] )
125 S_Y2 = variance( Y [start:stop+1] ) # to compute correlation
128 c = S_XY/(sqrt(S_X2)*sqrt(S_Y2))
129 a = S_XY/S_X2 # regr line coeffs
130 b= avg ( Y[start:stop+1] ) - a * avg( X[start:stop+1] )
131 print(" range [%d,%d] corr=%f, coeff det=%f [a=%f, b=%f]" % (X[start],X[stop],c,c**2,a, b))
132 correl.append( (c, stop-start) ); # store correl. coef + number of values (segment length)
133 interv.append( (a,b, X[start],X[stop]) );
136 glob_corr = glob_corr * c # product of correlation coeffs
137 print("-> glob_corr=%f\n" % glob_corr)
138 return (glob_corr,interv);
142 ##-----------------------------------------------------------------------------------------------
144 ##-----------------------------------------------------------------------------------------------
147 skampidat = open(sys.argv[1], "r")
151 for line in skampidat:
153 if line[0] != '#' and len(l)>=3: # is it a comment ?
157 #count= 8388608 8388608 144916.1 7.6 32 144916.1 143262.0
158 #("%s %d %d %f %f %d %f %f\n" % (countlbl, count, countn, time, stddev, iter, mini, maxi)
159 timings.append (float(l[3]))
160 sizes.append(int(l[1]))
165 ##----------------------- search for best break points-----------------
167 ## p1=2048 -> p1inx=11 delta=3 -> [8;14]
168 ## 8 : segments[(0,7),(8,13),(13,..)]
170 ## p2=65536 -> p2inx=16 delta=3 -> [13;19]
172 if len(sys.argv) == 4:
174 p1inx = sizes.index( p1 );
175 p2inx = sizes.index( p2 );
180 ## tweak parameters here to extend/reduce search
181 search_p1 = 30 # number of values to search +/- around p1
182 search_p2 = 45 # number of values to search +/- around p2
185 lb1 = max(1, p1inx-search_p1)
186 ub1 = min(p1inx+search_p1,search_p1, p2inx);
187 lb2 = max(p1inx,p2inx-search_p2) # breakpoint +/- delta
188 ub2 = min(p2inx+search_p2,len(sizes)-1);
190 print("** evaluating over \n");
191 print("interv1:\t %d <--- %d ---> %d" % (sizes[lb1],p1,sizes[ub1]))
192 print("rank: \t (%d)<---(%d)--->(%d)\n" % (lb1,p1inx,ub1))
193 print("interv2:\t\t %d <--- %d ---> %d" % (sizes[lb2],p2,sizes[ub2]))
194 print("rank: \t\t(%d)<---(%d)--->(%d)\n" % (lb2,p2inx,ub2))
195 for i in range(lb1,ub1+1):
196 for j in range(lb2,ub2+1):
197 if i<j: # segments must not overlap
198 if i+1 >=min_seg_size and j-i+1 >= min_seg_size and len(sizes)-1-j >= min_seg_size : # not too small segments
199 print("** i=%d,j=%d" % (i,j))
200 segments = [(0,i),(i,j),(j,len(sizes)-1)]
201 (glob_cor, interv) = correl_split( sizes, timings, segments)
202 if ( glob_cor > max_glob_corr):
203 max_glob_corr = glob_cor
206 for (a,b,i,j) in max_interv:
207 print("** OPT: [%d .. %d]" % (i,j))
208 print("** Product of correl coefs = %f" % (max_glob_corr))
210 print("#-------------------- cut here the gnuplot code -----------------------------------------------------------\n");
211 preamble='set output "regr.eps"\n\
212 set terminal postscript eps color\n\
214 set xlabel "Each message size in bytes"\n\
215 set ylabel "Time in us"\n\
221 print('plot "%s" u 3:4:($5) with errorbars title "skampi traces %s",\\' % (sys.argv[1],sys.argv[1]));
222 for (a,b,i,j) in max_interv:
223 print('"%s" u (%d<=$3 && $3<=%d? $3:0/0):(%f*($3)+%f) w linespoints title "regress. %s-%s bytes",\\' % (sys.argv[1],i,j,a,b,i,j))
225 print("#-------------------- /cut here the gnuplot code -----------------------------------------------------------\n");
229 print('\n** Linear regression on %d values **\n' % (nblines))
230 print('\n sizes=',sizes,'\n\n')
231 avg_sizes = avg( sizes )
232 avg_timings = avg( timings )
233 print("avg_timings=%f, avg_sizes=%f, nblines=%d\n" % (avg_timings,avg_sizes,nblines))
235 S_XY= cov( sizes, timings )
236 S_X2 = variance( sizes )
237 S_Y2 = variance( timings ) # to compute correlation
240 correl = S_XY/(sqrt(S_X2)*sqrt(S_Y2)) # corealation coeff (Bravais-Pearson)
243 b= avg_timings - a * avg_sizes
244 print("[S_XY=%f, S_X2=%f]\n[correlation=%f, coeff det=%f]\n[a=%f, b=%f]\n" % (S_XY, S_X2, correl,correl**2,a, b))