1 /* This file is part of the Linux Trace Toolkit viewer
2 * Copyright (C) 2009, 2010 Benjamin Poirier <benjamin.poirier@polymtl.ca>
4 * This program is free software: you can redistribute it and/or modify it
5 * under the terms of the GNU Lesser General Public License as published by
6 * the Free Software Foundation, either version 2.1 of the License, or (at
7 * your option) any later version.
9 * This program is distributed in the hope that it will be useful, but WITHOUT
10 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
12 * License for more details.
14 * You should have received a copy of the GNU Lesser General Public License
15 * along with this program. If not, see <http://www.gnu.org/licenses/>.
17 #define _ISOC99_SOURCE
32 #include "sync_chain.h"
34 #include "event_analysis_chull.h"
51 // Functions common to all analysis modules
52 static void initAnalysisCHull(SyncState
* const syncState
);
53 static void destroyAnalysisCHull(SyncState
* const syncState
);
55 static void analyzeMessageCHull(SyncState
* const syncState
, Message
* const
57 static AllFactors
* finalizeAnalysisCHull(SyncState
* const syncState
);
58 static void printAnalysisStatsCHull(SyncState
* const syncState
);
59 static void writeAnalysisGraphsPlotsCHull(SyncState
* const syncState
, const
60 unsigned int i
, const unsigned int j
);
62 // Functions specific to this module
63 static void openGraphFiles(SyncState
* const syncState
);
64 static void closeGraphFiles(SyncState
* const syncState
);
65 static void writeGraphFiles(SyncState
* const syncState
);
66 static void gfDumpHullToFile(gpointer data
, gpointer userData
);
68 static void grahamScan(GQueue
* const hull
, Point
* const newPoint
, const
70 static int jointCmp(const Point
* const p1
, const Point
* const p2
, const Point
*
71 const p3
) __attribute__((pure
));
72 static double crossProductK(const Point
const* p1
, const Point
const* p2
,
73 const Point
const* p3
, const Point
const* p4
) __attribute__((pure
));
74 static Factors
* calculateFactorsExact(GQueue
* const cu
, GQueue
* const cl
, const
75 LineType lineType
) __attribute__((pure
));
76 static void calculateFactorsFallback(GQueue
* const cr
, GQueue
* const cs
,
77 PairFactors
* const result
);
78 static double slope(const Point
* const p1
, const Point
* const p2
)
79 __attribute__((pure
));
80 static double intercept(const Point
* const p1
, const Point
* const p2
)
81 __attribute__((pure
));
82 static double verticalDistance(Point
* p1
, Point
* p2
, Point
* const point
)
83 __attribute__((pure
));
85 static void gfPointDestroy(gpointer data
, gpointer userData
);
88 static AnalysisModule analysisModuleCHull
= {
90 .initAnalysis
= &initAnalysisCHull
,
91 .destroyAnalysis
= &destroyAnalysisCHull
,
92 .analyzeMessage
= &analyzeMessageCHull
,
93 .finalizeAnalysis
= &finalizeAnalysisCHull
,
94 .printAnalysisStats
= &printAnalysisStatsCHull
,
96 .writeTraceTraceForePlots
= &writeAnalysisGraphsPlotsCHull
,
102 * Analysis module registering function
104 void registerAnalysisCHull()
106 g_queue_push_tail(&analysisModules
, &analysisModuleCHull
);
111 * Analysis init function
113 * This function is called at the beginning of a synchronization run for a set
116 * Allocate some of the analysis specific data structures
119 * syncState container for synchronization data.
120 * This function allocates or initializes these analysisData
125 static void initAnalysisCHull(SyncState
* const syncState
)
128 AnalysisDataCHull
* analysisData
;
130 analysisData
= malloc(sizeof(AnalysisDataCHull
));
131 syncState
->analysisData
= analysisData
;
133 analysisData
->hullArray
= malloc(syncState
->traceNb
* sizeof(GQueue
**));
134 for (i
= 0; i
< syncState
->traceNb
; i
++)
136 analysisData
->hullArray
[i
]= malloc(syncState
->traceNb
* sizeof(GQueue
*));
138 for (j
= 0; j
< syncState
->traceNb
; j
++)
140 analysisData
->hullArray
[i
][j
]= g_queue_new();
144 if (syncState
->stats
)
146 analysisData
->stats
= malloc(sizeof(AnalysisStatsCHull
));
147 analysisData
->stats
->dropped
= 0;
148 analysisData
->stats
->allFactors
= NULL
;
151 if (syncState
->graphsStream
)
153 analysisData
->graphsData
= malloc(sizeof(AnalysisGraphsDataCHull
));
154 openGraphFiles(syncState
);
155 analysisData
->graphsData
->allFactors
= NULL
;
161 * Create and open files used to store convex hull points to genereate
162 * graphs. Allocate and populate array to store file pointers.
165 * syncState: container for synchronization data
167 static void openGraphFiles(SyncState
* const syncState
)
173 AnalysisDataCHull
* analysisData
;
175 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
177 cwd
= changeToGraphsDir(syncState
->graphsDir
);
179 analysisData
->graphsData
->hullPoints
= malloc(syncState
->traceNb
*
181 for (i
= 0; i
< syncState
->traceNb
; i
++)
183 analysisData
->graphsData
->hullPoints
[i
]= malloc(syncState
->traceNb
*
185 for (j
= 0; j
< syncState
->traceNb
; j
++)
189 retval
= snprintf(name
, sizeof(name
),
190 "analysis_chull-%03u_to_%03u.data", j
, i
);
191 if (retval
> sizeof(name
) - 1)
193 name
[sizeof(name
) - 1]= '\0';
195 if ((analysisData
->graphsData
->hullPoints
[i
][j
]= fopen(name
, "w")) ==
198 g_error(strerror(errno
));
207 g_error(strerror(errno
));
214 * Write hull points to files to generate graphs.
217 * syncState: container for synchronization data
219 static void writeGraphFiles(SyncState
* const syncState
)
222 AnalysisDataCHull
* analysisData
;
224 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
226 for (i
= 0; i
< syncState
->traceNb
; i
++)
228 for (j
= 0; j
< syncState
->traceNb
; j
++)
232 g_queue_foreach(analysisData
->hullArray
[i
][j
],
234 analysisData
->graphsData
->hullPoints
[i
][j
]);
242 * A GFunc for g_queue_foreach. Write a hull point to a file used to generate
246 * data: Point*, point to write to the file
247 * userData: FILE*, file pointer where to write the point
249 static void gfDumpHullToFile(gpointer data
, gpointer userData
)
253 point
= (Point
*) data
;
254 fprintf((FILE*) userData
, "%20" PRIu64
" %20" PRIu64
"\n", point
->x
, point
->y
);
259 * Close files used to store convex hull points to generate graphs.
260 * Deallocate array to store file pointers.
263 * syncState: container for synchronization data
265 static void closeGraphFiles(SyncState
* const syncState
)
268 AnalysisDataCHull
* analysisData
;
271 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
273 if (analysisData
->graphsData
->hullPoints
== NULL
)
278 for (i
= 0; i
< syncState
->traceNb
; i
++)
280 for (j
= 0; j
< syncState
->traceNb
; j
++)
284 retval
= fclose(analysisData
->graphsData
->hullPoints
[i
][j
]);
287 g_error(strerror(errno
));
291 free(analysisData
->graphsData
->hullPoints
[i
]);
293 free(analysisData
->graphsData
->hullPoints
);
294 analysisData
->graphsData
->hullPoints
= NULL
;
299 * Analysis destroy function
301 * Free the analysis specific data structures
304 * syncState container for synchronization data.
305 * This function deallocates these analysisData members:
309 static void destroyAnalysisCHull(SyncState
* const syncState
)
312 AnalysisDataCHull
* analysisData
;
314 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
316 if (analysisData
== NULL
)
321 for (i
= 0; i
< syncState
->traceNb
; i
++)
323 for (j
= 0; j
< syncState
->traceNb
; j
++)
325 g_queue_foreach(analysisData
->hullArray
[i
][j
], gfPointDestroy
, NULL
);
326 g_queue_free(analysisData
->hullArray
[i
][j
]);
328 free(analysisData
->hullArray
[i
]);
330 free(analysisData
->hullArray
);
332 if (syncState
->stats
)
334 freeAllFactors(analysisData
->stats
->allFactors
);
336 free(analysisData
->stats
);
339 if (syncState
->graphsStream
)
341 if (analysisData
->graphsData
->hullPoints
!= NULL
)
343 closeGraphFiles(syncState
);
346 freeAllFactors(analysisData
->graphsData
->allFactors
);
348 free(analysisData
->graphsData
);
351 free(syncState
->analysisData
);
352 syncState
->analysisData
= NULL
;
357 * Perform analysis on an event pair.
360 * syncState container for synchronization data
361 * message structure containing the events
363 static void analyzeMessageCHull(SyncState
* const syncState
, Message
* const message
)
365 AnalysisDataCHull
* analysisData
;
370 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
372 newPoint
= malloc(sizeof(Point
));
373 if (message
->inE
->traceNum
< message
->outE
->traceNum
)
375 // CA is inE->traceNum
376 newPoint
->x
= message
->inE
->cpuTime
;
377 newPoint
->y
= message
->outE
->cpuTime
;
379 g_debug("Reception point hullArray[%lu][%lu] "
380 "x= inE->time= %" PRIu64
" y= outE->time= %" PRIu64
,
381 message
->inE
->traceNum
, message
->outE
->traceNum
, newPoint
->x
,
386 // CA is outE->traceNum
387 newPoint
->x
= message
->outE
->cpuTime
;
388 newPoint
->y
= message
->inE
->cpuTime
;
390 g_debug("Send point hullArray[%lu][%lu] "
391 "x= inE->time= %" PRIu64
" y= outE->time= %" PRIu64
,
392 message
->inE
->traceNum
, message
->outE
->traceNum
, newPoint
->x
,
397 analysisData
->hullArray
[message
->inE
->traceNum
][message
->outE
->traceNum
];
399 if (hull
->length
>= 1 && newPoint
->x
< ((Point
*)
400 g_queue_peek_tail(hull
))->x
)
402 if (syncState
->stats
)
404 analysisData
->stats
->dropped
++;
411 grahamScan(hull
, newPoint
, hullType
);
417 * Construct one half of a convex hull from abscissa-sorted points
420 * hull: the points already in the hull
421 * newPoint: a new point to consider
422 * type: which half of the hull to construct
424 static void grahamScan(GQueue
* const hull
, Point
* const newPoint
, const
429 g_debug("grahamScan(hull (length: %u), newPoint, %s)", hull
->length
, type
430 == LOWER
? "LOWER" : "UPPER");
441 if (hull
->length
>= 2)
443 g_debug("jointCmp(hull[%u], hull[%u], newPoint) * inversionFactor = %d * %d = %d",
446 jointCmp(g_queue_peek_nth(hull
, hull
->length
- 2),
447 g_queue_peek_tail(hull
), newPoint
),
449 jointCmp(g_queue_peek_nth(hull
, hull
->length
- 2),
450 g_queue_peek_tail(hull
), newPoint
) * inversionFactor
);
452 while (hull
->length
>= 2 && jointCmp(g_queue_peek_nth(hull
, hull
->length
-
453 2), g_queue_peek_tail(hull
), newPoint
) * inversionFactor
<= 0)
455 g_debug("Removing hull[%u]", hull
->length
);
456 free((Point
*) g_queue_pop_tail(hull
));
458 if (hull
->length
>= 2)
460 g_debug("jointCmp(hull[%u], hull[%u], newPoint) * inversionFactor = %d * %d = %d",
463 jointCmp(g_queue_peek_nth(hull
, hull
->length
- 2),
464 g_queue_peek_tail(hull
), newPoint
),
466 jointCmp(g_queue_peek_nth(hull
, hull
->length
- 2),
467 g_queue_peek_tail(hull
), newPoint
) * inversionFactor
);
470 g_queue_push_tail(hull
, newPoint
);
475 * Finalize the factor calculations
478 * syncState container for synchronization data.
481 * AllFactors* synchronization factors for each trace pair, the caller is
482 * responsible for freeing the structure
484 static AllFactors
* finalizeAnalysisCHull(SyncState
* const syncState
)
486 AnalysisDataCHull
* analysisData
;
487 AllFactors
* allFactors
;
489 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
491 if (syncState
->graphsStream
&& analysisData
->graphsData
->hullPoints
!= NULL
)
493 writeGraphFiles(syncState
);
494 closeGraphFiles(syncState
);
497 allFactors
= calculateAllFactors(syncState
);
499 if (syncState
->stats
)
501 allFactors
->refCount
++;
502 analysisData
->stats
->allFactors
= allFactors
;
505 if (syncState
->graphsStream
)
507 allFactors
->refCount
++;
508 analysisData
->graphsData
->allFactors
= allFactors
;
516 * Print statistics related to analysis. Must be called after
520 * syncState container for synchronization data.
522 static void printAnalysisStatsCHull(SyncState
* const syncState
)
524 AnalysisDataCHull
* analysisData
;
527 if (!syncState
->stats
)
532 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
534 printf("Convex hull analysis stats:\n");
535 printf("\tout of order packets dropped from analysis: %u\n",
536 analysisData
->stats
->dropped
);
538 printf("\tNumber of points in convex hulls:\n");
540 for (i
= 0; i
< syncState
->traceNb
; i
++)
542 for (j
= i
+ 1; j
< syncState
->traceNb
; j
++)
544 printf("\t\t%3d - %-3d: lower half-hull %-5u upper half-hull %-5u\n",
545 i
, j
, analysisData
->hullArray
[j
][i
]->length
,
546 analysisData
->hullArray
[i
][j
]->length
);
550 printf("\tIndividual synchronization factors:\n");
552 for (i
= 0; i
< syncState
->traceNb
; i
++)
554 for (j
= i
+ 1; j
< syncState
->traceNb
; j
++)
556 PairFactors
* factorsCHull
;
558 factorsCHull
= &analysisData
->stats
->allFactors
->pairFactors
[j
][i
];
559 printf("\t\t%3d - %-3d: %s", i
, j
,
560 approxNames
[factorsCHull
->type
]);
562 if (factorsCHull
->type
== EXACT
)
564 printf(" a0= % 7g a1= 1 %c %7g\n",
565 factorsCHull
->approx
->offset
,
566 factorsCHull
->approx
->drift
< 0. ? '-' : '+',
567 fabs(factorsCHull
->approx
->drift
));
569 else if (factorsCHull
->type
== ACCURATE
)
571 printf("\n\t\t a0: % 7g to % 7g (delta= %7g)\n",
572 factorsCHull
->max
->offset
, factorsCHull
->min
->offset
,
573 factorsCHull
->min
->offset
- factorsCHull
->max
->offset
);
574 printf("\t\t a1: 1 %+7g to %+7g (delta= %7g)\n",
575 factorsCHull
->min
->drift
- 1., factorsCHull
->max
->drift
-
576 1., factorsCHull
->max
->drift
- factorsCHull
->min
->drift
);
578 else if (factorsCHull
->type
== APPROXIMATE
)
580 printf(" a0= % 7g a1= 1 %c %7g error= %7g\n",
581 factorsCHull
->approx
->offset
, factorsCHull
->approx
->drift
582 - 1. < 0. ? '-' : '+', fabs(factorsCHull
->approx
->drift
-
583 1.), factorsCHull
->accuracy
);
585 else if (factorsCHull
->type
== INCOMPLETE
)
589 if (factorsCHull
->min
->drift
!= -INFINITY
)
591 printf("\t\t min: a0: % 7g a1: 1 %c %7g\n",
592 factorsCHull
->min
->offset
, factorsCHull
->min
->drift
-
593 1. < 0 ? '-' : '+', fabs(factorsCHull
->min
->drift
-
596 if (factorsCHull
->max
->drift
!= INFINITY
)
598 printf("\t\t max: a0: % 7g a1: 1 %c %7g\n",
599 factorsCHull
->max
->offset
, factorsCHull
->max
->drift
-
600 1. < 0 ? '-' : '+', fabs(factorsCHull
->max
->drift
-
604 else if (factorsCHull
->type
== SCREWED
)
608 if (factorsCHull
->min
!= NULL
&& factorsCHull
->min
->drift
!= -INFINITY
)
610 printf("\t\t min: a0: % 7g a1: 1 %c %7g\n",
611 factorsCHull
->min
->offset
, factorsCHull
->min
->drift
-
612 1. < 0 ? '-' : '+', fabs(factorsCHull
->min
->drift
-
615 if (factorsCHull
->max
!= NULL
&& factorsCHull
->max
->drift
!= INFINITY
)
617 printf("\t\t max: a0: % 7g a1: 1 %c %7g\n",
618 factorsCHull
->max
->offset
, factorsCHull
->max
->drift
-
619 1. < 0 ? '-' : '+', fabs(factorsCHull
->max
->drift
-
623 else if (factorsCHull
->type
== ABSENT
)
629 g_assert_not_reached();
637 * A GFunc for g_queue_foreach()
640 * data Point*, point to destroy
643 static void gfPointDestroy(gpointer data
, gpointer userData
)
647 point
= (Point
*) data
;
653 * Find out if a sequence of three points constitutes a "left turn" or a
657 * p1, p2, p3: The three points.
661 * 0 colinear (unlikely result since this uses floating point
665 static int jointCmp(const Point
const* p1
, const Point
const* p2
, const
669 const double fuzzFactor
= 0.;
671 result
= crossProductK(p1
, p2
, p1
, p3
);
672 g_debug("crossProductK(p1= (%" PRIu64
", %" PRIu64
"), "
673 "p2= (%" PRIu64
", %" PRIu64
"), p1= (%" PRIu64
", %" PRIu64
"), "
674 "p3= (%" PRIu64
", %" PRIu64
"))= %g",
675 p1
->x
, p1
->y
, p2
->x
, p2
->y
, p1
->x
, p1
->y
, p3
->x
, p3
->y
, result
);
676 if (result
< fuzzFactor
)
680 else if (result
> fuzzFactor
)
692 * Calculate the k component of the cross product of two vectors.
695 * p1, p2: start and end points of the first vector
696 * p3, p4: start and end points of the second vector
699 * the k component of the cross product when considering the two vectors to
700 * be in the i-j plane. The direction (sign) of the result can be useful to
701 * determine the relative orientation of the two vectors.
703 static double crossProductK(const Point
const* p1
, const Point
const* p2
,
704 const Point
const* p3
, const Point
const* p4
)
706 return ((double) p2
->x
- p1
->x
) * ((double) p4
->y
- p3
->y
) - ((double)
707 p2
->y
- p1
->y
) * ((double) p4
->x
- p3
->x
);
712 * Analyze the convex hulls to determine the synchronization factors between
713 * each pair of trace.
716 * syncState container for synchronization data.
719 * AllFactors*, see the documentation for the member allFactors of
720 * AnalysisStatsCHull.
722 AllFactors
* calculateAllFactors(SyncState
* const syncState
)
724 unsigned int traceNumA
, traceNumB
;
725 AllFactors
* allFactors
;
726 AnalysisDataCHull
* analysisData
;
728 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
730 // Allocate allFactors and calculate min and max
731 allFactors
= createAllFactors(syncState
->traceNb
);
732 for (traceNumA
= 0; traceNumA
< syncState
->traceNb
; traceNumA
++)
734 for (traceNumB
= 0; traceNumB
< traceNumA
; traceNumB
++)
741 size_t factorsOffset
;
743 {MINIMUM
, offsetof(PairFactors
, min
)},
744 {MAXIMUM
, offsetof(PairFactors
, max
)}
747 cr
= analysisData
->hullArray
[traceNumB
][traceNumA
];
748 cs
= analysisData
->hullArray
[traceNumA
][traceNumB
];
750 for (i
= 0; i
< sizeof(loopValues
) / sizeof(*loopValues
); i
++)
752 g_debug("allFactors[%u][%u].%s = calculateFactorsExact(cr= "
753 "hullArray[%u][%u], cs= hullArray[%u][%u], %s)",
754 traceNumA
, traceNumB
, loopValues
[i
].factorsOffset
==
755 offsetof(PairFactors
, min
) ? "min" : "max", traceNumB
,
756 traceNumA
, traceNumA
, traceNumB
, loopValues
[i
].lineType
==
757 MINIMUM
? "MINIMUM" : "MAXIMUM");
758 *((Factors
**) ((void*)
759 &allFactors
->pairFactors
[traceNumA
][traceNumB
] +
760 loopValues
[i
].factorsOffset
))=
761 calculateFactorsExact(cr
, cs
, loopValues
[i
].lineType
);
766 // Calculate approx when possible
767 for (traceNumA
= 0; traceNumA
< syncState
->traceNb
; traceNumA
++)
769 for (traceNumB
= 0; traceNumB
< traceNumA
; traceNumB
++)
771 PairFactors
* factorsCHull
;
773 factorsCHull
= &allFactors
->pairFactors
[traceNumA
][traceNumB
];
774 if (factorsCHull
->min
== NULL
&& factorsCHull
->max
== NULL
)
776 factorsCHull
->type
= APPROXIMATE
;
777 calculateFactorsFallback(analysisData
->hullArray
[traceNumB
][traceNumA
],
778 analysisData
->hullArray
[traceNumA
][traceNumB
],
779 &allFactors
->pairFactors
[traceNumA
][traceNumB
]);
781 else if (factorsCHull
->min
!= NULL
&& factorsCHull
->max
!= NULL
)
783 if (factorsCHull
->min
->drift
!= -INFINITY
&&
784 factorsCHull
->max
->drift
!= INFINITY
)
786 factorsCHull
->type
= ACCURATE
;
787 calculateFactorsMiddle(factorsCHull
);
789 else if (factorsCHull
->min
->drift
!= -INFINITY
||
790 factorsCHull
->max
->drift
!= INFINITY
)
792 factorsCHull
->type
= INCOMPLETE
;
796 factorsCHull
->type
= ABSENT
;
801 //g_assert_not_reached();
802 factorsCHull
->type
= SCREWED
;
811 /* Calculate approximative factors based on minimum and maximum limits. The
812 * best approximation to make is the interior bissector of the angle formed by
813 * the minimum and maximum lines.
815 * The formulae used come from [Haddad, Yoram: Performance dans les systèmes
816 * répartis: des outils pour les mesures, Université de Paris-Sud, Centre
817 * d'Orsay, September 1988] Section 6.1 p.44
819 * The reasoning for choosing this estimator comes from [Duda, A., Harrus, G.,
820 * Haddad, Y., and Bernard, G.: Estimating global time in distributed systems,
821 * Proc. 7th Int. Conf. on Distributed Computing Systems, Berlin, volume 18,
825 * factors: contains the min and max limits, used to store the result
827 void calculateFactorsMiddle(PairFactors
* const factors
)
829 double amin
, amax
, bmin
, bmax
, bhat
;
831 amin
= factors
->max
->offset
;
832 amax
= factors
->min
->offset
;
833 bmin
= factors
->min
->drift
;
834 bmax
= factors
->max
->drift
;
836 g_assert_cmpfloat(bmax
, >=, bmin
);
838 factors
->approx
= malloc(sizeof(Factors
));
839 bhat
= (bmax
* bmin
- 1. + sqrt(1. + pow(bmax
, 2.) * pow(bmin
, 2.) +
840 pow(bmax
, 2.) + pow(bmin
, 2.))) / (bmax
+ bmin
);
841 factors
->approx
->offset
= amax
- (amax
- amin
) / 2. * (pow(bhat
, 2.) + 1.)
842 / (1. + bhat
* bmax
);
843 factors
->approx
->drift
= bhat
;
844 factors
->accuracy
= bmax
- bmin
;
849 * Analyze the convex hulls to determine the minimum or maximum
850 * synchronization factors between one pair of trace.
852 * This implements and improves upon the algorithm in [Haddad, Yoram:
853 * Performance dans les systèmes répartis: des outils pour les mesures,
854 * Université de Paris-Sud, Centre d'Orsay, September 1988] Section 6.2 p.47
856 * Some degenerate cases are possible:
857 * 1) the result is unbounded. In that case, when searching for the maximum
858 * factors, result->drift= INFINITY; result->offset= -INFINITY. When
859 * searching for the minimum factors, it is the opposite. It is not
860 * possible to improve the situation with this data.
861 * 2) no line can be above the upper hull and below the lower hull. This is
862 * because the hulls intersect each other or are reversed. This means that
863 * an assertion was false. Most probably, the clocks are not linear. It is
864 * possible to repeat the search with another algorithm that will find a
865 * "best effort" approximation. See calculateFactorsApprox().
868 * cu: the upper half-convex hull, the line must pass above this
869 * and touch it in one point
870 * cl: the lower half-convex hull, the line must pass below this
871 * and touch it in one point
872 * lineType: search for minimum or maximum factors
875 * If a result is found, a struct Factors is allocated, filed with the
876 * result and returned
877 * NULL otherwise, degenerate case 2 is in effect
879 static Factors
* calculateFactorsExact(GQueue
* const cu
, GQueue
* const cl
, const
885 double inversionFactor
;
888 g_debug("calculateFactorsExact(cu= %p, cl= %p, %s)", cu
, cl
, lineType
==
889 MINIMUM
? "MINIMUM" : "MAXIMUM");
891 if (lineType
== MINIMUM
)
895 inversionFactor
= -1.;
907 // Check for degenerate case 1
908 if (c1
->length
== 0 || c2
->length
== 0 || ((Point
*) g_queue_peek_nth(c1
,
909 i1
))->x
>= ((Point
*) g_queue_peek_nth(c2
, i2
))->x
)
911 result
= malloc(sizeof(Factors
));
912 if (lineType
== MINIMUM
)
914 result
->drift
= -INFINITY
;
915 result
->offset
= INFINITY
;
919 result
->drift
= INFINITY
;
920 result
->offset
= -INFINITY
;
932 g_queue_peek_nth(c1
, i1
),
933 g_queue_peek_nth(c2
, i2
),
934 g_queue_peek_nth(c1
, i1
),
935 g_queue_peek_nth(c2
, i2
- 1)) * inversionFactor
< 0.
938 if (((Point
*) g_queue_peek_nth(c1
, i1
))->x
939 < ((Point
*) g_queue_peek_nth(c2
, i2
- 1))->x
)
951 i1
+ 1 < c1
->length
- 1
953 g_queue_peek_nth(c1
, i1
),
954 g_queue_peek_nth(c2
, i2
),
955 g_queue_peek_nth(c1
, i1
+ 1),
956 g_queue_peek_nth(c2
, i2
)) * inversionFactor
< 0.
959 if (((Point
*) g_queue_peek_nth(c1
, i1
+ 1))->x
960 < ((Point
*) g_queue_peek_nth(c2
, i2
))->x
)
974 g_queue_peek_nth(c1
, i1
),
975 g_queue_peek_nth(c2
, i2
),
976 g_queue_peek_nth(c1
, i1
),
977 g_queue_peek_nth(c2
, i2
- 1)) * inversionFactor
< 0.
980 p1
= g_queue_peek_nth(c1
, i1
);
981 p2
= g_queue_peek_nth(c2
, i2
);
983 g_debug("Resulting points are: c1[i1]: x= %" PRIu64
" y= %" PRIu64
984 " c2[i2]: x= %" PRIu64
" y= %" PRIu64
"", p1
->x
, p1
->y
, p2
->x
, p2
->y
);
986 result
= malloc(sizeof(Factors
));
987 result
->drift
= slope(p1
, p2
);
988 result
->offset
= intercept(p1
, p2
);
990 g_debug("Resulting factors are: drift= %g offset= %g", result
->drift
,
998 * Analyze the convex hulls to determine approximate synchronization factors
999 * between one pair of trace when there is no line that can fit in the
1000 * corridor separating them.
1002 * This implements the algorithm in [Ashton, P.: Algorithms for Off-line Clock
1003 * Synchronisation, University of Canterbury, December 1995, 26] Section 4.2.2
1006 * For each point p1 in cr
1007 * For each point p2 in cs
1009 * Calculate the line paramaters
1010 * For each point p3 in each convex hull
1011 * If p3 is on the wrong side of the line
1013 * If error < errorMin
1017 * cr: the upper half-convex hull
1018 * cs: the lower half-convex hull
1019 * result: a pointer to the pre-allocated struct where the results
1022 static void calculateFactorsFallback(GQueue
* const cr
, GQueue
* const cs
,
1023 PairFactors
* const result
)
1025 unsigned int i
, j
, k
;
1030 approx
= malloc(sizeof(Factors
));
1032 for (i
= 0; i
< cs
->length
; i
++)
1034 for (j
= 0; j
< cr
->length
; j
++)
1041 if (((Point
*) g_queue_peek_nth(cs
, i
))->x
< ((Point
*)g_queue_peek_nth(cr
, j
))->x
)
1043 p1
= *(Point
*)g_queue_peek_nth(cs
, i
);
1044 p2
= *(Point
*)g_queue_peek_nth(cr
, j
);
1048 p1
= *(Point
*)g_queue_peek_nth(cr
, j
);
1049 p2
= *(Point
*)g_queue_peek_nth(cs
, i
);
1052 // The lower hull should be above the point
1053 for (k
= 0; k
< cs
->length
; k
++)
1055 if (jointCmp(&p1
, &p2
, g_queue_peek_nth(cs
, k
)) < 0.)
1057 error
+= verticalDistance(&p1
, &p2
, g_queue_peek_nth(cs
, k
));
1061 // The upper hull should be below the point
1062 for (k
= 0; k
< cr
->length
; k
++)
1064 if (jointCmp(&p1
, &p2
, g_queue_peek_nth(cr
, k
)) > 0.)
1066 error
+= verticalDistance(&p1
, &p2
, g_queue_peek_nth(cr
, k
));
1070 if (error
< errorMin
)
1072 g_debug("Fallback: i= %u j= %u is a better match (error= %g)", i
, j
, error
);
1073 approx
->drift
= slope(&p1
, &p2
);
1074 approx
->offset
= intercept(&p1
, &p2
);
1080 result
->approx
= approx
;
1081 result
->accuracy
= errorMin
;
1086 * Calculate the vertical distance between a line and a point
1089 * p1, p2: Two points defining the line
1093 * the vertical distance
1095 static double verticalDistance(Point
* p1
, Point
* p2
, Point
* const point
)
1097 return fabs(slope(p1
, p2
) * point
->x
+ intercept(p1
, p2
) - point
->y
);
1102 * Calculate the slope between two points
1105 * p1, p2 the two points
1110 static double slope(const Point
* const p1
, const Point
* const p2
)
1112 return ((double) p2
->y
- p1
->y
) / (p2
->x
- p1
->x
);
1116 /* Calculate the y-intercept of a line that passes by two points
1119 * p1, p2 the two points
1124 static double intercept(const Point
* const p1
, const Point
* const p2
)
1126 return ((double) p2
->y
* p1
->x
- (double) p1
->y
* p2
->x
) / ((double) p1
->x
- p2
->x
);
1131 * Write the analysis-specific graph lines in the gnuplot script.
1134 * syncState: container for synchronization data
1135 * i: first trace number
1136 * j: second trace number, garanteed to be larger than i
1138 void writeAnalysisGraphsPlotsCHull(SyncState
* const syncState
, const unsigned
1139 int i
, const unsigned int j
)
1141 AnalysisDataCHull
* analysisData
;
1142 PairFactors
* factorsCHull
;
1144 analysisData
= (AnalysisDataCHull
*) syncState
->analysisData
;
1146 fprintf(syncState
->graphsStream
,
1147 "\t\"analysis_chull-%1$03d_to_%2$03d.data\" "
1148 "title \"Lower half-hull\" with linespoints "
1149 "linecolor rgb \"#015a01\" linetype 4 pointtype 8 pointsize 0.8, \\\n"
1150 "\t\"analysis_chull-%2$03d_to_%1$03d.data\" "
1151 "title \"Upper half-hull\" with linespoints "
1152 "linecolor rgb \"#003366\" linetype 4 pointtype 10 pointsize 0.8, \\\n",
1155 factorsCHull
= &analysisData
->graphsData
->allFactors
->pairFactors
[j
][i
];
1156 if (factorsCHull
->type
== EXACT
)
1158 fprintf(syncState
->graphsStream
,
1160 "title \"Exact conversion\" with lines "
1161 "linecolor rgb \"black\" linetype 1, \\\n",
1162 factorsCHull
->approx
->offset
, factorsCHull
->approx
->drift
);
1164 else if (factorsCHull
->type
== ACCURATE
)
1166 fprintf(syncState
->graphsStream
,
1167 "\t%.2f + %.10f * x "
1168 "title \"Min conversion\" with lines "
1169 "linecolor rgb \"black\" linetype 5, \\\n",
1170 factorsCHull
->min
->offset
, factorsCHull
->min
->drift
);
1171 fprintf(syncState
->graphsStream
,
1172 "\t%.2f + %.10f * x "
1173 "title \"Max conversion\" with lines "
1174 "linecolor rgb \"black\" linetype 8, \\\n",
1175 factorsCHull
->max
->offset
, factorsCHull
->max
->drift
);
1176 fprintf(syncState
->graphsStream
,
1177 "\t%.2f + %.10f * x "
1178 "title \"Middle conversion\" with lines "
1179 "linecolor rgb \"black\" linetype 1, \\\n",
1180 factorsCHull
->approx
->offset
, factorsCHull
->approx
->drift
);
1182 else if (factorsCHull
->type
== APPROXIMATE
)
1184 fprintf(syncState
->graphsStream
,
1185 "\t%.2f + %.10f * x "
1186 "title \"Fallback conversion\" with lines "
1187 "linecolor rgb \"gray60\" linetype 1, \\\n",
1188 factorsCHull
->approx
->offset
, factorsCHull
->approx
->drift
);
1190 else if (factorsCHull
->type
== INCOMPLETE
)
1192 if (factorsCHull
->min
->drift
!= -INFINITY
)
1194 fprintf(syncState
->graphsStream
,
1195 "\t%.2f + %.10f * x "
1196 "title \"Min conversion\" with lines "
1197 "linecolor rgb \"black\" linetype 5, \\\n",
1198 factorsCHull
->min
->offset
, factorsCHull
->min
->drift
);
1201 if (factorsCHull
->max
->drift
!= INFINITY
)
1203 fprintf(syncState
->graphsStream
,
1204 "\t%.2f + %.10f * x "
1205 "title \"Max conversion\" with lines "
1206 "linecolor rgb \"black\" linetype 8, \\\n",
1207 factorsCHull
->max
->offset
, factorsCHull
->max
->drift
);
1210 else if (factorsCHull
->type
== SCREWED
)
1212 if (factorsCHull
->min
!= NULL
&& factorsCHull
->min
->drift
!= -INFINITY
)
1214 fprintf(syncState
->graphsStream
,
1215 "\t%.2f + %.10f * x "
1216 "title \"Min conversion\" with lines "
1217 "linecolor rgb \"black\" linetype 5, \\\n",
1218 factorsCHull
->min
->offset
, factorsCHull
->min
->drift
);
1221 if (factorsCHull
->max
!= NULL
&& factorsCHull
->max
->drift
!= INFINITY
)
1223 fprintf(syncState
->graphsStream
,
1224 "\t%.2f + %.10f * x "
1225 "title \"Max conversion\" with lines "
1226 "linecolor rgb \"black\" linetype 8, \\\n",
1227 factorsCHull
->max
->offset
, factorsCHull
->max
->drift
);
1230 else if (factorsCHull
->type
== ABSENT
)
1235 g_assert_not_reached();