--- Begin Message ---
Package: release.debian.org
Severity: normal
User: release.debian.org@packages.debian.org
Usertags: unblock
Please unblock package ntop
Needed to remove non-DFSG code: http://bugs.debian.org/692732
Changes:
- replace CC-BY-NC countmin.[ch] and prng.[ch] with the new version
released under GPL-2. The difference is minimal, and mostly comments.
- the new countmin.[ch] and prng.[ch] are also included vanilla in the
debian directory so that get-orig-source.sh can create the new orig
tarball
- updated the Maintainer (moved myself from Uploader to Maintainer) as
agreed with former Maintainer
Debdiff attached.
unblock ntop/3:4.99.3+ndpi5517+dfsg2-1
-- System Information:
Debian Release: wheezy/sid
APT prefers testing
APT policy: (500, 'testing')
Architecture: amd64 (x86_64)
Foreign Architectures: i386
Kernel: Linux 3.2.0-4-amd64 (SMP w/2 CPU cores)
Locale: LANG=en_US.UTF-8, LC_CTYPE=en_US.UTF-8 (charmap=UTF-8)
Shell: /bin/sh linked to /bin/dash
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/countmin.c ntop-4.99.3+ndpi5517+dfsg2/countmin.c
--- ntop-4.99.3+ndpi5517+dfsg1/countmin.c 2011-11-21 05:42:08.000000000 -0800
+++ ntop-4.99.3+ndpi5517+dfsg2/countmin.c 2012-11-30 00:16:18.000000000 -0800
@@ -1,20 +1,30 @@
/********************************************************************
Count-Min Sketches
-G. Cormode 2003,2004
+G. Cormode 2003,2004, 2010, 2012
-Updated: 2004-06 Added a floating point sketch and support for
+Updated: 2004-06 Added a floating point sketch and support for
inner product point estimation
Initial version: 2003-12
-This work is licensed under the Creative Commons
-Attribution-NonCommercial License. To view a copy of this license,
-visit http://creativecommons.org/licenses/by-nc/1.0/ or send a letter
-to Creative Commons, 559 Nathan Abbott Way, Stanford, California
-94305, USA.
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License along
+ with this program; if not, write to the Free Software Foundation, Inc.,
+ 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
*********************************************************************/
#include <stdlib.h>
+#include "prng.h"
#include "countmin.h"
#define min(x,y) ((x) < (y) ? (x) : (y))
@@ -23,8 +33,6 @@
double eps; /* 1+epsilon = approximation factor */
double delta; /* probability of failure */
-//int bits=32;
-
/************************************************************************/
/* Routines to support Count-Min sketches */
/************************************************************************/
@@ -36,7 +44,7 @@
prng_type * prng;
cm=(CM_type *) malloc(sizeof(CM_type));
- prng=prng_Init(-abs(seed),2);
+ prng=prng_Init(-abs(seed),2);
// initialize the generator to pick the hash functions
if (cm && prng)
@@ -44,7 +52,6 @@
cm->depth=depth;
cm->width=width;
cm->count=0;
- cm->prng = prng; /* L. Deri */
cm->counts=(int **)calloc(sizeof(int *),cm->depth);
cm->counts[0]=(int *)calloc(sizeof(int), cm->depth*cm->width);
cm->hasha=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
@@ -61,6 +68,8 @@
}
else cm=NULL;
}
+ if (prng)
+ prng_Destroy(prng);
return cm;
}
@@ -105,7 +114,6 @@
}
if (cm->hasha) free(cm->hasha); cm->hasha=NULL;
if (cm->hashb) free(cm->hashb); cm->hashb=NULL;
- prng_Destroy(cm->prng); /* L.Deri */
free(cm); cm=NULL;
}
@@ -187,9 +195,10 @@
return 1;
}
-int CM_InnerProd(CM_type * cm1, CM_type * cm2)
+long long CM_InnerProd(CM_type * cm1, CM_type * cm2)
{ // Estimate the inner product of two vectors by comparing their sketches
- int i,j, tmp, result;
+ int i,j;
+ long long result, tmp;
result=0;
if (CM_Compatible(cm1,cm2))
@@ -207,9 +216,33 @@
return result;
}
+#if 0
+long long CM_F2Est(CM_type * cm)
+{ // Estimate the second frequency moment of the stream
+ int i,j;
+ long long result, tmp, *ans;
+
+ if (!cm) return 0;
+ ans=(long long *) calloc(1+cm->depth,sizeof(long long));
+
+ for (j=0;j<cm->depth;j++)
+ {
+ result=0;
+ for (i=0;i<cm->width;i+=2)
+ {
+ tmp=(cm->counts[j][i]-cm->counts[j][i+1]);
+ result+=tmp*tmp;
+ }
+ ans[j+1]=result;
+ }
+ result=LLMedSelect((cm->depth+1)/2,cm->depth,ans);
+ return result;
+}
+#endif
+
int CM_Residue(CM_type * cm, unsigned int * Q)
{
-// CM_Residue computes the sum of everything left after the points
+// CM_Residue computes the sum of everything left after the points
// from Q have been removed
// Q is a list of points, where Q[0] gives the length of the list
@@ -224,7 +257,7 @@
nextest=0;
for (i=0;i<cm->width;i++)
bitmap[i]=0;
- for (i=1;i<Q[0];i++)
+ for (i=1;i<(int) Q[0];i++)
bitmap[hash31(cm->hasha[j],cm->hashb[j],Q[i]) % cm->width]=1;
for (i=0;i<cm->width;i++)
if (bitmap[i]==0) nextest+=cm->counts[j][i];
@@ -245,7 +278,7 @@
cm=(CMF_type *) malloc(sizeof(CMF_type));
- prng=prng_Init(-abs(seed),2);
+ prng=prng_Init(-abs(seed),2);
// initialize the generator to pick the hash functions
if (cm && prng)
@@ -269,6 +302,8 @@
}
else cm=NULL;
}
+ if (prng)
+ prng_Destroy(prng);
return cm;
}
@@ -311,8 +346,6 @@
free(cm->counts);
cm->counts=NULL;
}
-
-
if (cm->hasha) free(cm->hasha); cm->hasha=NULL;
if (cm->hashb) free(cm->hashb); cm->hashb=NULL;
free(cm); cm=NULL;
@@ -339,10 +372,11 @@
// this can be done more efficiently if the width is a power of two
}
-int CMF_PointEst(CMF_type * cm, unsigned int query)
+double CMF_PointEst(CMF_type * cm, unsigned int query)
{
// return an estimate of the count of an item by taking the minimum
- int j, ans;
+ int j;
+ double ans=0.;
if (!cm) return 0;
ans=cm->counts[0][hash31(cm->hasha[0],cm->hashb[0],query) % cm->width];
@@ -380,12 +414,12 @@
{
loc=hash31(cm1->hasha[j],cm1->hashb[j],query) % cm1->width;
tmp=cm1->counts[j][loc]*cm2->counts[j][loc];
- ans=min(ans,tmp);
+ ans=min(ans,tmp);
}
}
return (ans);
}
-
+
double CMF_InnerProd(CMF_type * cm1, CMF_type * cm2)
{ // Estimate the inner product of two vectors by comparing their sketches
int i,j;
@@ -413,7 +447,7 @@
CMH_type * CMH_Init(int width, int depth, int U, int gran)
{
- // initialize a hierarchical set of sketches for range queries
+ // initialize a hierarchical set of sketches for range queries
// heavy hitters or quantiles
CMH_type * cmh;
@@ -424,10 +458,10 @@
// U is the log the size of the universe in bits
if (gran>U || gran<1) return(NULL);
- // gran is the granularity to look at the universe in
+ // gran is the granularity to look at the universe in
// check that the parameters make sense...
- cmh=(CMH_type *) malloc(sizeof(CMH_type));
+ cmh=(CMH_type *) calloc(1,sizeof(CMH_type));
prng=prng_Init(-12784,2);
// initialize the generator for picking the hash functions
@@ -441,11 +475,11 @@
cmh->gran=gran;
cmh->levels=(int) ceil(((float) U)/((float) gran));
for (j=0;j<cmh->levels;j++)
- if (1<<(cmh->gran*j) <= cmh->depth*cmh->width)
+ if ((long long) 1<<(cmh->gran*j) <= cmh->depth*cmh->width)
cmh->freelim=j;
//find the level up to which it is cheaper to keep exact counts
cmh->freelim=cmh->levels-cmh->freelim;
-
+
cmh->counts=(int **) calloc(sizeof(int *), 1+cmh->levels);
cmh->hasha=(unsigned int **)calloc(sizeof(unsigned int *),1+cmh->levels);
cmh->hashb=(unsigned int **)calloc(sizeof(unsigned int *),1+cmh->levels);
@@ -454,12 +488,12 @@
{
if (i>=cmh->freelim)
{ // allocate space for representing things exactly at high levels
- cmh->counts[i]=calloc(1<<(cmh->gran*j),sizeof(int));
+ cmh->counts[i]=(int *) calloc(1<<(cmh->gran*j),sizeof(int));
j++;
cmh->hasha[i]=NULL;
cmh->hashb[i]=NULL;
}
- else
+ else
{ // allocate space for a sketch
cmh->counts[i]=(int *)calloc(sizeof(int), cmh->depth*cmh->width);
cmh->hasha[i]=(unsigned int *)
@@ -476,11 +510,13 @@
}
}
}
+ if (prng)
+ prng_Destroy(prng);
return cmh;
}
void CMH_Destroy(CMH_type * cmh)
-{ // free up the space
+{ // free up the space
int i;
if (!cmh) return;
for (i=0;i<cmh->levels;i++)
@@ -489,7 +525,7 @@
{
free(cmh->counts[i]);
}
- else
+ else
{
free(cmh->hasha[i]);
free(cmh->hashb[i]);
@@ -513,14 +549,14 @@
{
offset=0;
if (i>=cmh->freelim)
- {
+ { // level 0 = leaves, higher levels = internal nodes
cmh->counts[i][item]+=diff;
- // keep exact counts at high levels in the hierarchy
+ // keep exact counts at high levels in the hierarchy
}
else
for (j=0;j<cmh->depth;j++)
{
- cmh->counts[i][(hash31(cmh->hasha[i][j],cmh->hashb[i][j],item)
+ cmh->counts[i][(hash31(cmh->hasha[i][j],cmh->hashb[i][j],item)
% cmh->width) + offset]+=diff;
// this can be done more efficiently if the width is a power of two
offset+=cmh->width;
@@ -545,7 +581,7 @@
return(admin + hashes + counts);
}
-int CMH_count(CMH_type * cmh, int depth, int item)
+int CMH_count(CMH_type * cmh, int depth, unsigned int item)
{
// return an estimate of item at level depth
@@ -561,23 +597,23 @@
// else, use the appropriate sketch to make an estimate
offset=0;
estimate=cmh->counts[depth][(hash31(cmh->hasha[depth][0],
- cmh->hashb[depth][0],item)
+ cmh->hashb[depth][0],item)
% cmh->width) + offset];
for (j=1;j<cmh->depth;j++)
{
offset+=cmh->width;
estimate=min(estimate,
cmh->counts[depth][(hash31(cmh->hasha[depth][j],
- cmh->hashb[depth][j],item)
+ cmh->hashb[depth][j],item)
% cmh->width) + offset]);
}
return(estimate);
}
-void CMH_recursive(CMH_type * cmh, int depth, int start,
+void CMH_recursive(CMH_type * cmh, int depth, int start,
int thresh, unsigned int * results)
{
- // for finding heavy hitters, recursively descend looking
+ // for finding heavy hitters, recursively descend looking
// for ranges that exceed the threshold
int i;
@@ -586,11 +622,11 @@
int itemshift;
estcount=CMH_count(cmh,depth,start);
- if (estcount>=thresh)
- {
+ if (estcount>=thresh)
+ {
if (depth==0)
{
- if (results[0]<cmh->width)
+ if (results[0]<(unsigned int) cmh->width)
{
results[0]++;
results[results[0]]=start;
@@ -607,7 +643,7 @@
}
}
-int * CMH_FindHH(CMH_type * cmh, int thresh)
+unsigned int * CMH_FindHH(CMH_type * cmh, int thresh)
{ // find all items whose estimated count is greater than phi n
unsigned int * results;
@@ -618,27 +654,31 @@
return(results);
}
-int CMH_Rangesum(CMH_type * cmh, int start, int end)
+int CMH_Rangesum(CMH_type * cmh, long long start, long long end)
{
- // compute a range sum:
+ // compute a range sum:
// start at bottom level
// compute any estimates needed at each level
// work upwards
- int leftend,rightend,i,depth, result, topend;
+ int depth, result;
+ long long range, leftend, rightend;
+ long long topend, i;
- topend=1<<cmh->U;
+ topend=((long long) 1)<<cmh->U;
end=min(topend,end);
- if ((end>topend) && (start==0))
+ if ((end>topend) && (start==0)) {
return cmh->count;
+ }
end+=1; // adjust for end effects
result=0;
for (depth=0;depth<=cmh->levels;depth++)
{
if (start==end) break;
- if ((end-start+1)<(1<<cmh->gran))
- { // at the highest level, avoid overcounting
+ range=(end-start+1);
+ if ((unsigned int) (end-start+1)<(((unsigned int) 1)<<cmh->gran))
+ { // if only a few nodes to probe at this level, probe them all
for (i=start;i<end;i++)
result+=CMH_count(cmh,depth,i);
break;
@@ -646,6 +686,7 @@
else
{ // figure out what needs to be done at each end
leftend=(((start>>cmh->gran)+1)<<cmh->gran) - start;
+ if (leftend>= 1<<cmh->gran) leftend=0;
rightend=(end)-((end>>cmh->gran)<<cmh->gran);
if ((leftend>0) && (start<end))
for (i=0;i<leftend;i++)
@@ -665,19 +706,20 @@
return result;
}
-int CMH_FindRange(CMH_type * cmh, int sum)
+long long CMH_FindRange(CMH_type * cmh, int sum)
{
- unsigned long low, high, mid=0, est;
- int i;
+ long long low, high, mid=0;
+ int i, est=0;
// find a range starting from zero that adds up to sum
if (cmh->count<sum) return 1<<(cmh->U);
low=0;
- high=1<<cmh->U;
+ high=((long long) 1)<<cmh->U;
+
for (i=0;i<cmh->U;i++)
{
mid=(low+high)/2;
- est=CMH_Rangesum(cmh,0,mid);
+ est=CMH_Rangesum(cmh,0,(unsigned int) mid);
if (est>sum)
high=mid;
else
@@ -687,18 +729,17 @@
}
-int CMH_AltFindRange(CMH_type * cmh, int sum)
+long long CMH_AltFindRange(CMH_type * cmh, int sum)
{
- unsigned long low, high, mid=0, est, top;
- int i;
+ long long low, high, mid=0,top;
+ int i, est=0;
// find a range starting from the right hand side that adds up to sum
if (cmh->count<sum) return 1<<(cmh->U);
low=0;
- top=1<<cmh->U;
+ top=((long long) 1)<<cmh->U;
high=top;
- for (i=0;i<cmh->U;i++)
- {
+ for (i=0;i<cmh->U;i++) {
mid=(low+high)/2;
est=CMH_Rangesum(cmh,mid,top);
if (est<sum)
@@ -710,17 +751,17 @@
}
-int CMH_Quantile(CMH_type * cmh, float frac)
+long long CMH_Quantile(CMH_type * cmh, float frac)
{
// find a quantile by doing the appropriate range search
if (frac<0) return 0;
- if (frac>1)
+ if (frac>1)
return 1<<cmh->U;
- return ((CMH_FindRange(cmh,cmh->count*frac)+
- CMH_AltFindRange(cmh,cmh->count*(1-frac)))/2);
+ return ((CMH_FindRange(cmh,(long long) (cmh->count*frac))+
+ CMH_AltFindRange(cmh,(long long) (cmh->count*(1-frac))))/2);
// each result gives a lower/upper bound on the location of the quantile
// with high probability, these will be close: only a small number of values
- // will be between the estimates.
+ // will be between the estimates.
}
long long CMH_F2Est(CMH_type * cmh)
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/countmin.h ntop-4.99.3+ndpi5517+dfsg2/countmin.h
--- ntop-4.99.3+ndpi5517+dfsg1/countmin.h 2012-02-21 15:51:02.000000000 -0800
+++ ntop-4.99.3+ndpi5517+dfsg2/countmin.h 2012-11-30 00:16:18.000000000 -0800
@@ -2,15 +2,8 @@
// 1 -- The basic CM Sketch
// 2 -- The hierarchical CM Sketch: with log n levels, for range sums etc.
-#ifndef min
#define min(x,y) ((x) < (y) ? (x) : (y))
-#endif
-
-#ifndef max
#define max(x,y) ((x) > (y) ? (x) : (y))
-#endif
-
-#include "prng.h"
typedef struct CM_type{
long long count;
@@ -18,7 +11,6 @@
int width;
int ** counts;
unsigned int *hasha, *hashb;
- prng_type * prng; /* L. Deri */
} CM_type;
typedef struct CMF_type{ // shadow of above stucture with floats
@@ -37,8 +29,9 @@
extern void CM_Update(CM_type *, unsigned int, int);
extern int CM_PointEst(CM_type *, unsigned int);
extern int CM_PointMed(CM_type *, unsigned int);
-extern int CM_InnerProd(CM_type *, CM_type *);
+extern long long CM_InnerProd(CM_type *, CM_type *);
extern int CM_Residue(CM_type *, unsigned int *);
+extern long long CM_F2Est(CM_type *);
extern CMF_type * CMF_Init(int, int, int);
extern CMF_type * CMF_Copy(CMF_type *);
@@ -57,7 +50,7 @@
int depth;
int width;
int ** counts;
- unsigned int **hasha, * *hashb;
+ unsigned int **hasha, **hashb;
} CMH_type;
extern CMH_type * CMH_Init(int, int, int, int);
@@ -66,9 +59,9 @@
extern int CMH_Size(CMH_type *);
extern void CMH_Update(CMH_type *, unsigned int, int);
-extern int * CMH_FindHH(CMH_type *, int);
-extern int CMH_Rangesum(CMH_type *, int, int);
+extern unsigned int * CMH_FindHH(CMH_type *, int);
+extern int CMH_Rangesum(CMH_type *, long long, long long);
-extern int CMH_FindRange(CMH_type * cmh, int);
-extern int CMH_Quantile(CMH_type *cmh,float);
+extern long long CMH_FindRange(CMH_type * cmh, int);
+extern long long CMH_Quantile(CMH_type *cmh,float);
extern long long CMH_F2Est(CMH_type *);
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/debian/changelog ntop-4.99.3+ndpi5517+dfsg2/debian/changelog
--- ntop-4.99.3+ndpi5517+dfsg1/debian/changelog 2012-06-27 00:25:38.000000000 -0700
+++ ntop-4.99.3+ndpi5517+dfsg2/debian/changelog 2012-11-30 00:37:15.000000000 -0800
@@ -1,3 +1,11 @@
+ntop (3:4.99.3+ndpi5517+dfsg2-1) unstable; urgency=medium
+
+ * Repackage upstream source replacing non-DFSG countmin code with the GPL
+ one (Closes: #692732).
+ * Update Maintainer and Uploaders.
+
+ -- Ludovico Cavedon <cavedon@debian.org> Fri, 30 Nov 2012 00:07:10 -0800
+
ntop (3:4.99.3+ndpi5517+dfsg1-1) unstable; urgency=low
* Imported Upstream version 4.99.3 and nDPI r5517.
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/debian/control ntop-4.99.3+ndpi5517+dfsg2/debian/control
--- ntop-4.99.3+ndpi5517+dfsg1/debian/control 2012-06-27 00:25:38.000000000 -0700
+++ ntop-4.99.3+ndpi5517+dfsg2/debian/control 2012-11-30 00:37:15.000000000 -0800
@@ -1,8 +1,7 @@
Source: ntop
Section: net
Priority: optional
-Maintainer: Jordan Metzmeier <jordan@linuxgen.com>
-Uploaders: Ludovico Cavedon <cavedon@debian.org>
+Maintainer: Ludovico Cavedon <cavedon@debian.org>
Build-Depends-Indep: groff
Build-Depends: debhelper (>= 8), libpcap0.8-dev, autoconf, automake,
libgdbm-dev, libpng-dev, libncurses5-dev, libtool, libwrap0-dev,
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/debian/copyright ntop-4.99.3+ndpi5517+dfsg2/debian/copyright
--- ntop-4.99.3+ndpi5517+dfsg1/debian/copyright 2012-06-27 00:25:38.000000000 -0700
+++ ntop-4.99.3+ndpi5517+dfsg2/debian/copyright 2012-11-30 00:37:15.000000000 -0800
@@ -57,13 +57,9 @@
WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED WARRANTIES OF
MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
-Files: countmin.h
-Copyright: 2003-2004, G. Cormode
-License: CC-BY-NC
- This work is licensed under the Creative Commons Attribution-NonCommercial
- License. To view a copy of this license, visit
- http://creativecommons.org/licenses/by-nc/1.0/ or send a letter to Creative
- Commons, 559 Nathan Abbott Way, Stanford, California 94305, USA.
+Files: countmin.h countmin.c prng.h prng.c
+Copyright: 2003-2004, 2010, 2012, G. Cormode
+License: GPL-2+
Files: utils/xmldump.awk utils/processstruct.awk AutoConfigureNtop linuxrelease
config_h.shm utildl.c
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/debian/countmin.c ntop-4.99.3+ndpi5517+dfsg2/debian/countmin.c
--- ntop-4.99.3+ndpi5517+dfsg1/debian/countmin.c 1969-12-31 16:00:00.000000000 -0800
+++ ntop-4.99.3+ndpi5517+dfsg2/debian/countmin.c 2012-11-30 00:37:15.000000000 -0800
@@ -0,0 +1,788 @@
+/********************************************************************
+Count-Min Sketches
+
+G. Cormode 2003,2004, 2010, 2012
+
+Updated: 2004-06 Added a floating point sketch and support for
+ inner product point estimation
+Initial version: 2003-12
+
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License along
+ with this program; if not, write to the Free Software Foundation, Inc.,
+ 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+*********************************************************************/
+
+#include <stdlib.h>
+#include "prng.h"
+#include "countmin.h"
+
+#define min(x,y) ((x) < (y) ? (x) : (y))
+#define max(x,y) ((x) > (y) ? (x) : (y))
+
+double eps; /* 1+epsilon = approximation factor */
+double delta; /* probability of failure */
+
+/************************************************************************/
+/* Routines to support Count-Min sketches */
+/************************************************************************/
+
+CM_type * CM_Init(int width, int depth, int seed)
+{ // Initialize the sketch based on user-supplied size
+ CM_type * cm;
+ int j;
+ prng_type * prng;
+
+ cm=(CM_type *) malloc(sizeof(CM_type));
+ prng=prng_Init(-abs(seed),2);
+ // initialize the generator to pick the hash functions
+
+ if (cm && prng)
+ {
+ cm->depth=depth;
+ cm->width=width;
+ cm->count=0;
+ cm->counts=(int **)calloc(sizeof(int *),cm->depth);
+ cm->counts[0]=(int *)calloc(sizeof(int), cm->depth*cm->width);
+ cm->hasha=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
+ cm->hashb=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
+ if (cm->counts && cm->hasha && cm->hashb && cm->counts[0])
+ {
+ for (j=0;j<depth;j++)
+ {
+ cm->hasha[j]=prng_int(prng) & MOD;
+ cm->hashb[j]=prng_int(prng) & MOD;
+ // pick the hash functions
+ cm->counts[j]=(int *) cm->counts[0]+(j*cm->width);
+ }
+ }
+ else cm=NULL;
+ }
+ if (prng)
+ prng_Destroy(prng);
+ return cm;
+}
+
+CM_type * CM_Copy(CM_type * cmold)
+{ // create a new sketch with the same parameters as an existing one
+ CM_type * cm;
+ int j;
+
+ if (!cmold) return(NULL);
+ cm=(CM_type *) malloc(sizeof(CM_type));
+ if (cm)
+ {
+ cm->depth=cmold->depth;
+ cm->width=cmold->width;
+ cm->count=0;
+ cm->counts=(int **)calloc(sizeof(int *),cm->depth);
+ cm->counts[0]=(int *)calloc(sizeof(int), cm->depth*cm->width);
+ cm->hasha=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
+ cm->hashb=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
+ if (cm->counts && cm->hasha && cm->hashb && cm->counts[0])
+ {
+ for (j=0;j<cm->depth;j++)
+ {
+ cm->hasha[j]=cmold->hasha[j];
+ cm->hashb[j]=cmold->hashb[j];
+ cm->counts[j]=(int *) cm->counts[0]+(j*cm->width);
+ }
+ }
+ else cm=NULL;
+ }
+ return cm;
+}
+
+void CM_Destroy(CM_type * cm)
+{ // get rid of a sketch and free up the space
+ if (!cm) return;
+ if (cm->counts)
+ {
+ if (cm->counts[0]) free(cm->counts[0]);
+ free(cm->counts);
+ cm->counts=NULL;
+ }
+ if (cm->hasha) free(cm->hasha); cm->hasha=NULL;
+ if (cm->hashb) free(cm->hashb); cm->hashb=NULL;
+ free(cm); cm=NULL;
+}
+
+int CM_Size(CM_type * cm)
+{ // return the size of the sketch in bytes
+ int counts, hashes, admin;
+ if (!cm) return 0;
+ admin=sizeof(CM_type);
+ counts=cm->width*cm->depth*sizeof(int);
+ hashes=cm->depth*2*sizeof(unsigned int);
+ return(admin + hashes + counts);
+}
+
+void CM_Update(CM_type * cm, unsigned int item, int diff)
+{
+ int j;
+
+ if (!cm) return;
+ cm->count+=diff;
+ for (j=0;j<cm->depth;j++)
+ cm->counts[j][hash31(cm->hasha[j],cm->hashb[j],item) % cm->width]+=diff;
+ // this can be done more efficiently if the width is a power of two
+}
+
+int CM_PointEst(CM_type * cm, unsigned int query)
+{
+ // return an estimate of the count of an item by taking the minimum
+ int j, ans;
+
+ if (!cm) return 0;
+ ans=cm->counts[0][hash31(cm->hasha[0],cm->hashb[0],query) % cm->width];
+ for (j=1;j<cm->depth;j++)
+ ans=min(ans,cm->counts[j][hash31(cm->hasha[j],cm->hashb[j],query)%cm->width]);
+ // this can be done more efficiently if the width is a power of two
+ return (ans);
+}
+
+#if 0
+int CM_PointMed(CM_type * cm, unsigned int query)
+{
+ // return an estimate of the count by taking the median estimate
+ // useful when counts can become negative
+ // depth needs to be larger for this to work well
+ int j, * ans, result=0;
+
+ if (!cm) return 0;
+ ans=(int *) calloc(1+cm->depth,sizeof(int));
+ for (j=0;j<cm->depth;j++)
+ ans[j+1]=cm->counts[j][hash31(cm->hasha[j],cm->hashb[j],query)%cm->width];
+
+ if (cm->depth==1)
+ result=ans[1];
+ else
+ if (cm->depth==2)
+ {
+ //result=(ans[1]+ans[2])/2;
+ if (abs(ans[1]) < abs(ans[2]))
+ result=ans[1]; else result=ans[2];
+ // special tweak for small depth sketches
+ }
+ else
+ result=(MedSelect(1+cm->depth/2,cm->depth,ans));
+ return result;
+ // need to adjust for routine starting at 1
+}
+#endif
+
+int CM_Compatible(CM_type * cm1, CM_type * cm2)
+{ // test whether two sketches are comparable (have same parameters)
+ int i;
+ if (!cm1 || !cm2) return 0;
+ if (cm1->width!=cm2->width) return 0;
+ if (cm1->depth!=cm2->depth) return 0;
+ for (i=0;i<cm1->depth;i++)
+ {
+ if (cm1->hasha[i]!=cm2->hasha[i]) return 0;
+ if (cm1->hashb[i]!=cm2->hashb[i]) return 0;
+ }
+ return 1;
+}
+
+long long CM_InnerProd(CM_type * cm1, CM_type * cm2)
+{ // Estimate the inner product of two vectors by comparing their sketches
+ int i,j;
+ long long result, tmp;
+
+ result=0;
+ if (CM_Compatible(cm1,cm2))
+ {
+ for (i=0;i<cm1->width;i++)
+ result+=cm1->counts[0][i]*cm2->counts[0][i];
+ for (j=1;j<cm1->depth;j++)
+ {
+ tmp=0;
+ for (i=0;i<cm1->width;i++)
+ tmp+=cm1->counts[j][i]*cm2->counts[j][i];
+ result=min(tmp,result);
+ }
+ }
+ return result;
+}
+
+#if 0
+long long CM_F2Est(CM_type * cm)
+{ // Estimate the second frequency moment of the stream
+ int i,j;
+ long long result, tmp, *ans;
+
+ if (!cm) return 0;
+ ans=(long long *) calloc(1+cm->depth,sizeof(long long));
+
+ for (j=0;j<cm->depth;j++)
+ {
+ result=0;
+ for (i=0;i<cm->width;i+=2)
+ {
+ tmp=(cm->counts[j][i]-cm->counts[j][i+1]);
+ result+=tmp*tmp;
+ }
+ ans[j+1]=result;
+ }
+ result=LLMedSelect((cm->depth+1)/2,cm->depth,ans);
+ return result;
+}
+#endif
+
+int CM_Residue(CM_type * cm, unsigned int * Q)
+{
+// CM_Residue computes the sum of everything left after the points
+// from Q have been removed
+// Q is a list of points, where Q[0] gives the length of the list
+
+ char * bitmap;
+ int i,j;
+ int estimate=0, nextest;
+
+ if (!cm) return 0;
+ bitmap=(char *) calloc(cm->width,sizeof(char));
+ for (j=0;j<cm->depth;j++)
+ {
+ nextest=0;
+ for (i=0;i<cm->width;i++)
+ bitmap[i]=0;
+ for (i=1;i<(int) Q[0];i++)
+ bitmap[hash31(cm->hasha[j],cm->hashb[j],Q[i]) % cm->width]=1;
+ for (i=0;i<cm->width;i++)
+ if (bitmap[i]==0) nextest+=cm->counts[j][i];
+ estimate=max(estimate,nextest);
+ }
+ return(estimate);
+}
+
+/************************************************************************/
+/* Routines to support Count-Min sketches with floating point data */
+/************************************************************************/
+
+CMF_type * CMF_Init(int width, int depth, int seed)
+{ // Initialize the sketch based on user-supplied size
+ CMF_type * cm;
+ int j;
+ prng_type * prng;
+
+ cm=(CMF_type *) malloc(sizeof(CMF_type));
+
+ prng=prng_Init(-abs(seed),2);
+ // initialize the generator to pick the hash functions
+
+ if (cm && prng)
+ {
+ cm->depth=depth;
+ cm->width=width;
+ cm->count=0;
+ cm->counts=(double **)calloc(sizeof(double *),cm->depth);
+ cm->counts[0]=(double *)calloc(sizeof(double), cm->depth*cm->width);
+ cm->hasha=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
+ cm->hashb=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
+ if (cm->counts && cm->hasha && cm->hashb && cm->counts[0])
+ {
+ for (j=0;j<depth;j++)
+ {
+ cm->hasha[j]=prng_int(prng) & MOD;
+ cm->hashb[j]=prng_int(prng) & MOD;
+ // pick the hash functions
+ cm->counts[j]=(double *) cm->counts[0]+(j*cm->width);
+ }
+ }
+ else cm=NULL;
+ }
+ if (prng)
+ prng_Destroy(prng);
+ return cm;
+}
+
+CMF_type * CMF_Copy(CMF_type * cmold)
+{ // create a new sketch with the same parameters as an existing one
+ CMF_type * cm;
+ int j;
+
+ if (!cmold) return(NULL);
+ cm=(CMF_type *) malloc(sizeof(CMF_type));
+ if (cm)
+ {
+ cm->depth=cmold->depth;
+ cm->width=cmold->width;
+ cm->count=0;
+ cm->counts=(double **)calloc(sizeof(double *),cm->depth);
+ cm->counts[0]=(double *)calloc(sizeof(double), cm->depth*cm->width);
+ cm->hasha=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
+ cm->hashb=(unsigned int *)calloc(sizeof(unsigned int),cm->depth);
+ if (cm->counts && cm->hasha && cm->hashb && cm->counts[0])
+ {
+ for (j=0;j<cm->depth;j++)
+ {
+ cm->hasha[j]=cmold->hasha[j];
+ cm->hashb[j]=cmold->hashb[j];
+ cm->counts[j]=(double *) cm->counts[0]+(j*cm->width);
+ }
+ }
+ else cm=NULL;
+ }
+ return cm;
+}
+
+void CMF_Destroy(CMF_type * cm)
+{ // get rid of a sketch and free up the space
+ if (!cm) return;
+ if (cm->counts)
+ {
+ if (cm->counts[0]) free(cm->counts[0]);
+ free(cm->counts);
+ cm->counts=NULL;
+ }
+ if (cm->hasha) free(cm->hasha); cm->hasha=NULL;
+ if (cm->hashb) free(cm->hashb); cm->hashb=NULL;
+ free(cm); cm=NULL;
+}
+
+int CMF_Size(CMF_type * cm)
+{ // return the size of the sketch in bytes
+ int counts, hashes, admin;
+ if (!cm) return 0;
+ admin=sizeof(CM_type);
+ counts=cm->width*cm->depth*sizeof(double);
+ hashes=cm->depth*2*sizeof(unsigned int);
+ return(admin + hashes + counts);
+}
+
+void CMF_Update(CMF_type * cm, unsigned int item, double diff)
+{
+ int j;
+
+ if (!cm) return;
+ cm->count+=diff;
+ for (j=0;j<cm->depth;j++)
+ cm->counts[j][hash31(cm->hasha[j],cm->hashb[j],item) % cm->width]+=diff;
+ // this can be done more efficiently if the width is a power of two
+}
+
+double CMF_PointEst(CMF_type * cm, unsigned int query)
+{
+ // return an estimate of the count of an item by taking the minimum
+ int j;
+ double ans=0.;
+
+ if (!cm) return 0;
+ ans=cm->counts[0][hash31(cm->hasha[0],cm->hashb[0],query) % cm->width];
+ for (j=1;j<cm->depth;j++)
+ ans=min(ans,cm->counts[j][hash31(cm->hasha[j],cm->hashb[j],query)%cm->width]);
+ // this can be done more efficiently if the width is a power of two
+ return (ans);
+}
+
+int CMF_Compatible(CMF_type * cm1, CMF_type * cm2)
+{ // test whether two sketches are comparable (have same parameters)
+ int i;
+ if (!cm1 || !cm2) return 0;
+ if (cm1->width!=cm2->width) return 0;
+ if (cm1->depth!=cm2->depth) return 0;
+ for (i=0;i<cm1->depth;i++)
+ {
+ if (cm1->hasha[i]!=cm2->hasha[i]) return 0;
+ if (cm1->hashb[i]!=cm2->hashb[i]) return 0;
+ }
+ return 1;
+}
+
+double CMF_PointProd(CMF_type * cm1, CMF_type * cm2, unsigned int query)
+{ // Estimate the inner product of two vectors by comparing their sketches
+ int j, loc;
+ double tmp, ans;
+
+ ans=0.0;
+ if (CMF_Compatible(cm1,cm2))
+ {
+ loc=hash31(cm1->hasha[0],cm1->hashb[0],query) % cm1->width;
+ ans=cm1->counts[0][loc]*cm2->counts[0][loc];
+ for (j=1;j<cm1->depth;j++)
+ {
+ loc=hash31(cm1->hasha[j],cm1->hashb[j],query) % cm1->width;
+ tmp=cm1->counts[j][loc]*cm2->counts[j][loc];
+ ans=min(ans,tmp);
+ }
+ }
+ return (ans);
+}
+
+double CMF_InnerProd(CMF_type * cm1, CMF_type * cm2)
+{ // Estimate the inner product of two vectors by comparing their sketches
+ int i,j;
+ double tmp, result;
+
+ result=0;
+ if (CMF_Compatible(cm1,cm2))
+ {
+ for (i=0;i<cm1->width;i++)
+ result+=cm1->counts[0][i]*cm2->counts[0][i];
+ for (j=1;j<cm1->depth;j++)
+ {
+ tmp=0.0;
+ for (i=0;i<cm1->width;i++)
+ tmp+=cm1->counts[j][i]*cm2->counts[j][i];
+ result=min(tmp,result);
+ }
+ }
+ return result;
+}
+
+/************************************************************************/
+/* Routines to support hierarchical Count-Min sketches */
+/************************************************************************/
+
+CMH_type * CMH_Init(int width, int depth, int U, int gran)
+{
+ // initialize a hierarchical set of sketches for range queries
+ // heavy hitters or quantiles
+
+ CMH_type * cmh;
+ int i,j,k;
+ prng_type * prng;
+
+ if (U<=0 || U>32) return(NULL);
+ // U is the log the size of the universe in bits
+
+ if (gran>U || gran<1) return(NULL);
+ // gran is the granularity to look at the universe in
+ // check that the parameters make sense...
+
+ cmh=(CMH_type *) calloc(1,sizeof(CMH_type));
+
+ prng=prng_Init(-12784,2);
+ // initialize the generator for picking the hash functions
+
+ if (cmh && prng)
+ {
+ cmh->depth=depth;
+ cmh->width=width;
+ cmh->count=0;
+ cmh->U=U;
+ cmh->gran=gran;
+ cmh->levels=(int) ceil(((float) U)/((float) gran));
+ for (j=0;j<cmh->levels;j++)
+ if ((long long) 1<<(cmh->gran*j) <= cmh->depth*cmh->width)
+ cmh->freelim=j;
+ //find the level up to which it is cheaper to keep exact counts
+ cmh->freelim=cmh->levels-cmh->freelim;
+
+ cmh->counts=(int **) calloc(sizeof(int *), 1+cmh->levels);
+ cmh->hasha=(unsigned int **)calloc(sizeof(unsigned int *),1+cmh->levels);
+ cmh->hashb=(unsigned int **)calloc(sizeof(unsigned int *),1+cmh->levels);
+ j=1;
+ for (i=cmh->levels-1;i>=0;i--)
+ {
+ if (i>=cmh->freelim)
+ { // allocate space for representing things exactly at high levels
+ cmh->counts[i]=(int *) calloc(1<<(cmh->gran*j),sizeof(int));
+ j++;
+ cmh->hasha[i]=NULL;
+ cmh->hashb[i]=NULL;
+ }
+ else
+ { // allocate space for a sketch
+ cmh->counts[i]=(int *)calloc(sizeof(int), cmh->depth*cmh->width);
+ cmh->hasha[i]=(unsigned int *)
+ calloc(sizeof(unsigned int),cmh->depth);
+ cmh->hashb[i]=(unsigned int *)
+ calloc(sizeof(unsigned int),cmh->depth);
+
+ if (cmh->hasha[i] && cmh->hashb[i])
+ for (k=0;k<cmh->depth;k++)
+ { // pick the hash functions
+ cmh->hasha[i][k]=prng_int(prng) & MOD;
+ cmh->hashb[i][k]=prng_int(prng) & MOD;
+ }
+ }
+ }
+ }
+ if (prng)
+ prng_Destroy(prng);
+ return cmh;
+}
+
+void CMH_Destroy(CMH_type * cmh)
+{ // free up the space
+ int i;
+ if (!cmh) return;
+ for (i=0;i<cmh->levels;i++)
+ {
+ if (i>=cmh->freelim)
+ {
+ free(cmh->counts[i]);
+ }
+ else
+ {
+ free(cmh->hasha[i]);
+ free(cmh->hashb[i]);
+ free(cmh->counts[i]);
+ }
+ }
+ free(cmh->counts);
+ free(cmh->hasha);
+ free(cmh->hashb);
+ free(cmh);
+ cmh=NULL;
+}
+
+void CMH_Update(CMH_type * cmh, unsigned int item, int diff)
+{ // update with a new value
+ int i,j,offset;
+
+ if (!cmh) return;
+ cmh->count+=diff;
+ for (i=0;i<cmh->levels;i++)
+ {
+ offset=0;
+ if (i>=cmh->freelim)
+ { // level 0 = leaves, higher levels = internal nodes
+ cmh->counts[i][item]+=diff;
+ // keep exact counts at high levels in the hierarchy
+ }
+ else
+ for (j=0;j<cmh->depth;j++)
+ {
+ cmh->counts[i][(hash31(cmh->hasha[i][j],cmh->hashb[i][j],item)
+ % cmh->width) + offset]+=diff;
+ // this can be done more efficiently if the width is a power of two
+ offset+=cmh->width;
+ }
+ item>>=cmh->gran;
+ }
+}
+
+int CMH_Size(CMH_type * cmh)
+{ // return the size used in bytes
+ int counts, hashes, admin,i;
+ if (!cmh) return 0;
+ admin=sizeof(CMH_type);
+ counts=cmh->levels*sizeof(int **);
+ for (i=0;i<cmh->levels;i++)
+ if (i>=cmh->freelim)
+ counts+=(1<<(cmh->gran*(cmh->levels-i)))*sizeof(int);
+ else
+ counts+=cmh->width*cmh->depth*sizeof(int);
+ hashes=(cmh->levels-cmh->freelim)*cmh->depth*2*sizeof(unsigned int);
+ hashes+=(cmh->levels)*sizeof(unsigned int *);
+ return(admin + hashes + counts);
+}
+
+int CMH_count(CMH_type * cmh, int depth, unsigned int item)
+{
+ // return an estimate of item at level depth
+
+ int j;
+ int offset;
+ int estimate;
+
+ if (depth>=cmh->levels) return(cmh->count);
+ if (depth>=cmh->freelim)
+ { // use an exact count if there is one
+ return(cmh->counts[depth][item]);
+ }
+ // else, use the appropriate sketch to make an estimate
+ offset=0;
+ estimate=cmh->counts[depth][(hash31(cmh->hasha[depth][0],
+ cmh->hashb[depth][0],item)
+ % cmh->width) + offset];
+ for (j=1;j<cmh->depth;j++)
+ {
+ offset+=cmh->width;
+ estimate=min(estimate,
+ cmh->counts[depth][(hash31(cmh->hasha[depth][j],
+ cmh->hashb[depth][j],item)
+ % cmh->width) + offset]);
+ }
+ return(estimate);
+}
+
+void CMH_recursive(CMH_type * cmh, int depth, int start,
+ int thresh, unsigned int * results)
+{
+ // for finding heavy hitters, recursively descend looking
+ // for ranges that exceed the threshold
+
+ int i;
+ int blocksize;
+ int estcount;
+ int itemshift;
+
+ estcount=CMH_count(cmh,depth,start);
+ if (estcount>=thresh)
+ {
+ if (depth==0)
+ {
+ if (results[0]<(unsigned int) cmh->width)
+ {
+ results[0]++;
+ results[results[0]]=start;
+ }
+ }
+ else
+ {
+ blocksize=1<<cmh->gran;
+ itemshift=start<<cmh->gran;
+ // assumes that gran is an exact multiple of the bit dept
+ for (i=0;i<blocksize;i++)
+ CMH_recursive(cmh,depth-1,itemshift+i,thresh,results);
+ }
+ }
+}
+
+unsigned int * CMH_FindHH(CMH_type * cmh, int thresh)
+{ // find all items whose estimated count is greater than phi n
+
+ unsigned int * results;
+ results=(unsigned int *) calloc(cmh->width,sizeof(unsigned int));
+ results[0]=0;
+
+ CMH_recursive(cmh,cmh->levels,0,thresh,results);
+ return(results);
+}
+
+int CMH_Rangesum(CMH_type * cmh, long long start, long long end)
+{
+ // compute a range sum:
+ // start at bottom level
+ // compute any estimates needed at each level
+ // work upwards
+
+ int depth, result;
+ long long range, leftend, rightend;
+ long long topend, i;
+
+ topend=((long long) 1)<<cmh->U;
+ end=min(topend,end);
+ if ((end>topend) && (start==0)) {
+ return cmh->count;
+ }
+
+ end+=1; // adjust for end effects
+ result=0;
+ for (depth=0;depth<=cmh->levels;depth++)
+ {
+ if (start==end) break;
+ range=(end-start+1);
+ if ((unsigned int) (end-start+1)<(((unsigned int) 1)<<cmh->gran))
+ { // if only a few nodes to probe at this level, probe them all
+ for (i=start;i<end;i++)
+ result+=CMH_count(cmh,depth,i);
+ break;
+ }
+ else
+ { // figure out what needs to be done at each end
+ leftend=(((start>>cmh->gran)+1)<<cmh->gran) - start;
+ if (leftend>= 1<<cmh->gran) leftend=0;
+ rightend=(end)-((end>>cmh->gran)<<cmh->gran);
+ if ((leftend>0) && (start<end))
+ for (i=0;i<leftend;i++)
+ {
+ result+=CMH_count(cmh,depth,start+i);
+ }
+ if ((rightend>0) && (start<end))
+ for (i=0;i<rightend;i++)
+ {
+ result+=CMH_count(cmh,depth,end-i-1);
+ }
+ start=start>>cmh->gran;
+ if (leftend>0) start++;
+ end=end>>cmh->gran;
+ }
+ }
+ return result;
+}
+
+long long CMH_FindRange(CMH_type * cmh, int sum)
+{
+ long long low, high, mid=0;
+ int i, est=0;
+ // find a range starting from zero that adds up to sum
+
+ if (cmh->count<sum) return 1<<(cmh->U);
+ low=0;
+ high=((long long) 1)<<cmh->U;
+
+ for (i=0;i<cmh->U;i++)
+ {
+ mid=(low+high)/2;
+ est=CMH_Rangesum(cmh,0,(unsigned int) mid);
+ if (est>sum)
+ high=mid;
+ else
+ low=mid;
+ }
+ return mid;
+
+}
+
+long long CMH_AltFindRange(CMH_type * cmh, int sum)
+{
+ long long low, high, mid=0,top;
+ int i, est=0;
+ // find a range starting from the right hand side that adds up to sum
+
+ if (cmh->count<sum) return 1<<(cmh->U);
+ low=0;
+ top=((long long) 1)<<cmh->U;
+ high=top;
+ for (i=0;i<cmh->U;i++) {
+ mid=(low+high)/2;
+ est=CMH_Rangesum(cmh,mid,top);
+ if (est<sum)
+ high=mid;
+ else
+ low=mid;
+ }
+ return mid;
+
+}
+
+long long CMH_Quantile(CMH_type * cmh, float frac)
+{
+ // find a quantile by doing the appropriate range search
+ if (frac<0) return 0;
+ if (frac>1)
+ return 1<<cmh->U;
+ return ((CMH_FindRange(cmh,(long long) (cmh->count*frac))+
+ CMH_AltFindRange(cmh,(long long) (cmh->count*(1-frac))))/2);
+ // each result gives a lower/upper bound on the location of the quantile
+ // with high probability, these will be close: only a small number of values
+ // will be between the estimates.
+}
+
+long long CMH_F2Est(CMH_type * cmh)
+{
+ // A heuristic for estimating the F2 of a stream
+ // tends to overestimate a great deal on non-skewed streams
+
+ int i,j,k;
+ long long est, result;
+
+ k=0; result=-1;
+ for (i=0;i<cmh->depth;i++)
+ {
+ est=0;
+ for (j=0;j<cmh->width;j++)
+ {
+ est+=(long long) cmh->counts[0][k] * (long long) cmh->counts[0][k];
+ k++;
+ }
+ if (result<0) result=est; else
+ result=min(result,est);
+ }
+ return result;
+}
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/debian/countmin.h ntop-4.99.3+ndpi5517+dfsg2/debian/countmin.h
--- ntop-4.99.3+ndpi5517+dfsg1/debian/countmin.h 1969-12-31 16:00:00.000000000 -0800
+++ ntop-4.99.3+ndpi5517+dfsg2/debian/countmin.h 2012-11-30 00:37:15.000000000 -0800
@@ -0,0 +1,67 @@
+// Two different structures:
+// 1 -- The basic CM Sketch
+// 2 -- The hierarchical CM Sketch: with log n levels, for range sums etc.
+
+#define min(x,y) ((x) < (y) ? (x) : (y))
+#define max(x,y) ((x) > (y) ? (x) : (y))
+
+typedef struct CM_type{
+ long long count;
+ int depth;
+ int width;
+ int ** counts;
+ unsigned int *hasha, *hashb;
+} CM_type;
+
+typedef struct CMF_type{ // shadow of above stucture with floats
+ double count;
+ int depth;
+ int width;
+ double ** counts;
+ unsigned int *hasha, *hashb;
+} CMF_type;
+
+extern CM_type * CM_Init(int, int, int);
+extern CM_type * CM_Copy(CM_type *);
+extern void CM_Destroy(CM_type *);
+extern int CM_Size(CM_type *);
+
+extern void CM_Update(CM_type *, unsigned int, int);
+extern int CM_PointEst(CM_type *, unsigned int);
+extern int CM_PointMed(CM_type *, unsigned int);
+extern long long CM_InnerProd(CM_type *, CM_type *);
+extern int CM_Residue(CM_type *, unsigned int *);
+extern long long CM_F2Est(CM_type *);
+
+extern CMF_type * CMF_Init(int, int, int);
+extern CMF_type * CMF_Copy(CMF_type *);
+extern void CMF_Destroy(CMF_type *);
+extern int CMF_Size(CMF_type *);
+extern void CMF_Update(CMF_type *, unsigned int, double);
+extern double CMF_InnerProd(CMF_type *, CMF_type *);
+extern double CMF_PointProd(CMF_type *, CMF_type *, unsigned int);
+
+typedef struct CMH_type{
+ long long count;
+ int U; // size of the universe in bits
+ int gran; // granularity: eg 1, 4 or 8 bits
+ int levels; // function of U and gran
+ int freelim; // up to which level to keep exact counts
+ int depth;
+ int width;
+ int ** counts;
+ unsigned int **hasha, **hashb;
+} CMH_type;
+
+extern CMH_type * CMH_Init(int, int, int, int);
+extern CMH_type * CMH_Copy(CMH_type *);
+extern void CMH_Destroy(CMH_type *);
+extern int CMH_Size(CMH_type *);
+
+extern void CMH_Update(CMH_type *, unsigned int, int);
+extern unsigned int * CMH_FindHH(CMH_type *, int);
+extern int CMH_Rangesum(CMH_type *, long long, long long);
+
+extern long long CMH_FindRange(CMH_type * cmh, int);
+extern long long CMH_Quantile(CMH_type *cmh,float);
+extern long long CMH_F2Est(CMH_type *);
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/debian/get-orig-source.sh ntop-4.99.3+ndpi5517+dfsg2/debian/get-orig-source.sh
--- ntop-4.99.3+ndpi5517+dfsg1/debian/get-orig-source.sh 2012-06-27 00:25:38.000000000 -0700
+++ ntop-4.99.3+ndpi5517+dfsg2/debian/get-orig-source.sh 2012-11-30 00:37:15.000000000 -0800
@@ -8,6 +8,8 @@
UPSTREAM_DIR=$(echo -n $DEB_UPSTREAM_VERSION | sed -e 's/+rc[0-9]*//')
+curdir=$(dirname $0)
+
if [ ! -f ../ntop_$DEB_SOURCE_VERSION.orig.tar.gz ]; then
uscan --noconf --force-download --rename --download-current-version --destdir=.
rm -rf ntop-$DEB_UPSTREAM_VERSION
@@ -23,6 +25,13 @@
# remove binary-onyl files from MaxMind
rm ntop-$UPSTREAM_DIR/3rd_party/Geo*.dat.gz
+ # replace non-DFSG countmin with the GPL one
+ rm ntop-$UPSTREAM_DIR/countmin.[ch]
+ rm ntop-$UPSTREAM_DIR/prng.[ch]
+
+ cp $curdir/countmin.[ch] ntop-$UPSTREAM_DIR/
+ cp $curdir/prng.[ch] ntop-$UPSTREAM_DIR/
+
# remove non-DFSG-compliant part of ntop.h
ed ntop-$UPSTREAM_DIR/ntop.h > /dev/null <<EOF
/Declaration of POSIX directory browsing functions and types for Win32.
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/debian/prng.c ntop-4.99.3+ndpi5517+dfsg2/debian/prng.c
--- ntop-4.99.3+ndpi5517+dfsg1/debian/prng.c 1969-12-31 16:00:00.000000000 -0800
+++ ntop-4.99.3+ndpi5517+dfsg2/debian/prng.c 2012-11-30 00:37:15.000000000 -0800
@@ -0,0 +1,525 @@
+/*********************************************************************
+
+prng: pseudo-random number generators and hashing routines
+
+G. Cormode 2003-2012
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License along
+ with this program; if not, write to the Free Software Foundation, Inc.,
+ 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+************************************************************************/
+
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include "prng.h"
+
+#define PI 3.141592653589793
+
+long hash31(long long a, long long b, long long x)
+{
+
+ long long result;
+ long lresult;
+
+ // return a hash of x using a and b mod (2^31 - 1)
+// may need to do another mod afterwards, or drop high bits
+// depending on d, number of bad guys
+// 2^31 - 1 = 2147483647
+
+ // result = ((long long) a)*((long long) x)+((long long) b);
+ result=(a * x) + b;
+ result = ((result >> HL) + result) & MOD;
+ lresult=(long) result;
+
+ return(lresult);
+}
+
+long fourwise(long long a, long long b, long long c, long long d, long long x)
+{
+ long long result;
+ long lresult;
+
+ // returns values that are 4-wise independent by repeated calls
+ // to the pairwise indpendent routine.
+
+ result = hash31(hash31(hash31(a,b,x),c,x),d,x);
+ lresult = (long) result;
+ return lresult;
+}
+
+
+/*************************************************************************/
+/* First, some pseudo-random number generators sourced from other places */
+/*************************************************************************/
+
+// There are *THREE* alternate implementations of PRNGs here.
+// One taken from Numerical Recipes in C, the second from www.agner.org
+// The third is an internal C random library, srand
+// The variable usenric controls which one is used: pick one
+// and stick with it, switching between the two will give unpredictable
+// results. This is controlled by the randinit procedure, call it with
+// usenric == 1 to use the Numerical Recipes gens
+// usenric == 2 to use the agner.org PRNGs or
+// usenric == 3 to use the inbuilt C routines
+
+// from the math library:
+extern double sqrt(double);
+
+// following definitions needed for the random number generator
+#define IA 16807
+#define IM 2147483647
+#define AM (1.0/IM)
+#define IQ 127773
+#define IR 2836
+#define NDIV (1+(IM-1)/NTAB)
+#define EPS 1.2e-7
+#define RNMX (1.0-EPS)
+
+float ran1(prng_type * prng) {
+
+ // A Random Number Generator that picks a uniform [0,1] random number
+ // From Numerical Recipes, page 280
+ // Should be called with a NEGATIVE value of idum to initialize
+ // subsequent calls should not alter idum
+
+ int j;
+ long k;
+ float temp;
+
+ if (prng->floatidum <= 0 || !prng->iy) {
+ if (-(prng->floatidum) < 1) prng->floatidum=1;
+ else prng->floatidum = -(prng->floatidum);
+ for (j=NTAB+7;j>=0;j--) {
+ k=(prng->floatidum)/IQ;
+ prng->floatidum=IA*(prng->floatidum-k*IQ)-IR*k;
+ if (prng->floatidum < 0) prng->floatidum+=IM;
+ if (j<NTAB) prng->iv[j]=prng->floatidum;
+ }
+ prng->iy=prng->iv[0];
+ }
+ k = (prng->floatidum)/IQ;
+ prng->floatidum=IA*(prng->floatidum-k*IQ)-IR*k;
+ if (prng->floatidum<0) prng->floatidum += IM;
+ j = prng->iy/NDIV;
+ prng->iy=prng->iv[j];
+ prng->iv[j]=prng->floatidum;
+ if ((temp=AM*prng->iy) > RNMX) return RNMX;
+ else return temp;
+}
+
+long ran2(prng_type * prng) {
+
+ // A Random Number Generator that picks a uniform random number
+ // from the range of long integers.
+ // From Numerical Recipes, page 280
+ // Should be called with a NEGATIVE value of idum to initialize
+ // subsequent calls should not alter idum
+ // This is a hacked version of the above procedure, so proceed with
+ // caution.
+
+ int j;
+ long k;
+
+ if (prng->intidum <= 0 || !prng->iy) {
+ if (-(prng->intidum) < 1) prng->intidum=1;
+ else prng->intidum = -(prng->intidum);
+ for (j=NTAB+7;j>=0;j--) {
+ k=(prng->intidum)/IQ;
+ prng->intidum=IA*(prng->intidum-k*IQ)-IR*k;
+ if (prng->intidum < 0) prng->intidum+=IM;
+ if (j<NTAB) prng->iv[j]=prng->intidum;
+ }
+ prng->iy=prng->iv[0];
+ }
+ k = (prng->intidum)/IQ;
+ prng->intidum=IA*(prng->intidum-k*IQ)-IR*k;
+ if (prng->intidum<0) prng->intidum += IM;
+ j = prng->iy/NDIV;
+ prng->iy=prng->iv[j];
+ prng->iv[j]=prng->intidum;
+ return prng->iy;
+}
+
+/**********************************************************************/
+
+// Following routines are from www.agner.org
+
+/************************* RANROTB.C ******************** AgF 1999-03-03 *
+* Random Number generator 'RANROT' type B *
+* *
+* This is a lagged-Fibonacci type of random number generator with *
+* rotation of bits. The algorithm is: *
+* X[n] = ((X[n-j] rotl r1) + (X[n-k] rotl r2)) modulo 2^b *
+* *
+* The last k values of X are stored in a circular buffer named *
+* randbuffer. *
+* *
+* This version works with any integer size: 16, 32, 64 bits etc. *
+* The integers must be unsigned. The resolution depends on the integer *
+* size. *
+* *
+* Note that the function RanrotAInit must be called before the first *
+* call to RanrotA or iRanrotA *
+* *
+* The theory of the RANROT type of generators is described at *
+* www.agner.org/random/ranrot.htm *
+* *
+*************************************************************************/
+
+// this should be almost verbatim from the above webpage.
+// although it's been hacked with a little bit...
+
+unsigned long rotl (unsigned long x, unsigned long r) {
+ return (x << r) | (x >> (sizeof(x)*8-r));}
+
+/* define parameters (R1 and R2 must be smaller than the integer size): */
+#define JJ 10
+#define R1 5
+#define R2 3
+
+/* returns some random bits */
+unsigned long ran3(prng_type * prng) {
+
+ unsigned long x;
+
+ /* generate next random number */
+
+ x = prng->randbuffer[prng->r_p1] = rotl(prng->randbuffer[prng->r_p2], R1)
+ + rotl(prng->randbuffer[prng->r_p1], R2);
+ /* rotate list pointers */
+ if (--prng->r_p1 < 0) prng->r_p1 = KK - 1;
+ if (--prng->r_p2 < 0) prng->r_p2 = KK - 1;
+ /* conversion to float */
+ return x;
+}
+
+/* returns a random number between 0 and 1 */
+double ran4(prng_type * prng) {
+
+ /* conversion to floating point type */
+ return (ran3(prng) * prng->scale);
+}
+
+/* this function initializes the random number generator. */
+/* Must be called before the first call to RanrotA or iRanrotA */
+void RanrotAInit (prng_type * prng, unsigned long seed) {
+
+ int i;
+
+ /* put semi-random numbers into the buffer */
+ for (i=0; i<KK; i++) {
+ prng->randbuffer[i] = seed;
+ seed = rotl(seed,5) + 97;}
+
+ /* initialize pointers to circular buffer */
+ prng->r_p1 = 0; prng->r_p2 = JJ;
+
+ /* randomize */
+ for (i = 0; i < 300; i++) ran3(prng);
+ prng->scale = ldexp(1, -8*sizeof(unsigned long));
+}
+
+
+/**********************************************************************/
+/* These are wrapper procedures for the uniform random number gens */
+/**********************************************************************/
+
+long prng_int(prng_type * prng) {
+
+ // returns a pseudo-random long integer. Initialise the generator
+ // before use!
+
+ long response=0;
+
+ switch (prng->usenric)
+ {
+ case 1 : response=(ran2(prng)); break;
+ case 2 : response=(ran3(prng)); break;
+ case 3 : response=(lrand48()); break;
+ }
+ return response;
+}
+
+
+float prng_float(prng_type * prng) {
+
+ // returns a pseudo-random float in the range [0.0,1.0].
+ // Initialise the generator before use!
+ float result=0;
+
+ switch (prng->usenric)
+ {
+ case 1 : result=(ran1(prng)); break;
+ case 2 : result=(ran4(prng)); break;
+ case 3 : result=(drand48()); break;
+ }
+ return result;
+}
+
+prng_type * prng_Init(long seed, int nric) {
+
+ // Initialise the random number generators. nric determines
+ // which algorithm to use, 1 for Numerical Recipes in C,
+ // 0 for the other one.
+ prng_type * result;
+
+ result=(prng_type *) calloc(1,sizeof(prng_type));
+
+ result->iy=0;
+ result->usenric=nric;
+ result->floatidum=-1;
+ result->intidum=-1;
+ result->iset=0;
+ // set a global variable to record which algorithm to use
+ switch (nric)
+ {
+ case 2 :
+ RanrotAInit(result,seed);
+ break;
+ case 1 :
+ if (seed>0) {
+ // to initialise the NRiC PRNGs, call it with a negative value
+ // so make sure it gets a negative value!
+ result->floatidum = -(seed); result->intidum = -(seed);
+ } else {
+ result->floatidum=seed; result->intidum=seed;
+ }
+ break;
+ case 3 :
+ srand48(seed);
+ break;
+ }
+
+ prng_float(result);
+ prng_int(result);
+ // call the routines to actually initialise them
+ return(result);
+}
+
+void prng_Reseed(prng_type * prng, long seed)
+{
+ switch (prng->usenric)
+ {
+ case 2 :
+ RanrotAInit(prng,seed);
+ break;
+ case 1 :
+ if (seed>0) {
+ // to initialise the NRiC PRNGs, call it with a negative value
+ // so make sure it gets a negative value!
+ prng->floatidum = -(seed); prng->intidum = -(seed);
+ } else {
+ prng->floatidum=seed; prng->intidum=seed;
+ }
+ break;
+ case 3 :
+ srand48(seed);
+ break;
+ }
+}
+
+void prng_Destroy(prng_type * prng)
+{
+ free(prng);
+}
+
+/**********************************************************************/
+/* Next, a load of routines that convert uniform random variables */
+/* from [0,1] to stable distribitions, such as gaussian, levy or */
+/* general */
+/**********************************************************************/
+
+double prng_normal(prng_type * prng) {
+
+ // Pick random values distributed N(0,1) using the Box-Muller transform
+ // Taken from numerical recipes in C p289
+ // picks two at a time, returns one per call (buffers the other)
+
+ double fac,rsq,v1,v2;
+ if (prng->iset == 0) {
+ do {
+ v1 = 2.0*prng_float(prng)-1.0;
+ v2 = 2.0*prng_float(prng)-1.0;
+ rsq=v1*v1+v2*v2;
+ } while (rsq >= 1.0 || rsq == 0.0);
+
+ fac = sqrt((double) -2.0*log((double)rsq)/rsq);
+ prng->gset=v1*fac;
+ prng->iset=1;
+ return v2*fac;
+ }
+ else {
+ prng->iset = 0;
+ return prng->gset;
+ }
+}
+
+double prng_stabledbn(prng_type * prng, double alpha) {
+
+ // From 'stable distributions', John Nolan, manuscript, p24
+ // we set beta = 0 by analogy with the normal and cauchy case
+ // identical to the above routine, but returns a double instead
+ // of a long double (you'll see this a lot...)
+
+ double theta, W, holder, left, right;
+
+ theta=PI*(prng_float(prng) - 0.5);
+ W = -log(prng_float(prng)); // takes natural log
+
+ // printf("theta %Lf, W = %Lf \n", theta, W);
+
+ // some notes on Nolan's notes:
+ // if beta == 0 then c(alpha,beta)=1; theta_0 = 0
+ // expression reduces to sin alpha.theta / (cos theta) ^1/alpha
+ // * (cos (theta - alpha theta)/W) ^(1-alpha)/alpha
+ left = (sin(alpha*theta)/pow(cos(theta), 1.0/alpha));
+ right= pow(cos(theta*(1.0 - alpha))/W, ((1.0-alpha)/alpha));
+ holder=left*right;
+ return(holder);
+}
+
+
+long double prng_cauchy(prng_type * prng) {
+
+ // return a value from the cauchy distribution
+ // using the formula in 'Stable Distributions', p23
+ // this is distributed Cauchy(1,0)
+
+ return(tan(PI*(prng_float(prng) - 0.5)));
+}
+
+
+double prng_altstab(prng_type * prng, double p)
+{
+ double u,v,result;
+
+ u=prng_float(prng);
+ v=prng_float(prng);
+ result=pow(u,p);
+ // result=exp(p*log(u));
+ if (v<0.5) result=-result;
+ return(result);
+}
+
+
+/*
+long double levy() {
+
+ // this would give the levy distribution, except it doesn't get used
+
+ long double z;
+
+ z=gasdev();
+ return (1.0/(z*z));
+ }
+
+*/
+
+double prng_stable(prng_type * prng, double alpha) {
+
+ // wrapper for the stable distributions above:
+ // call the appropriate routine based on the value of alpha given
+ // initialising it with the seed in idum
+
+ // randinit must be called before entering this procedure for
+ // the first time since it uses the random generators
+
+
+ if (alpha==2.0)
+ return(prng_normal(prng));
+ else if (alpha==1.0)
+ return(prng_cauchy(prng));
+ else if (alpha<0.01)
+ return(prng_altstab(prng,-50.0));
+ else return (prng_stabledbn(prng,alpha));
+}
+
+double zeta(long n, double theta)
+{
+
+ // the zeta function, used by the below zipf function
+ // (this is not often called from outside this library)
+ // ... but have made it public now to speed things up
+
+ int i;
+ double ans=0.0;
+
+ for (i=1; i <= n; i++)
+ ans += pow(1./(double)i, theta);
+ return(ans);
+}
+
+double fastzipf(double theta, long n, double zetan, prng_type * prng) {
+
+ // this draws values from the zipf distribution
+ // this is mainly useful for test generation purposes
+ // n is range, theta is skewness parameter
+ // theta = 0 gives uniform dbn,
+ // theta > 1 gives highly skewed dbn.
+ // original code due to Flip Korn, used with permission
+
+ double alpha;
+ double eta;
+ double u;
+ double uz;
+ long double val;
+
+ // randinit must be called before entering this procedure for
+ // the first time since it uses the random generators
+
+ alpha = 1. / (1. - theta);
+ eta = (1. - pow(2./((double) n), 1. - theta))
+ / (1. - zeta(2,theta)/zetan);
+
+ u = prng_float(prng);
+ uz = u * zetan;
+ if (uz < 1.) val = 1;
+ else if (uz < (1. + pow(0.5, theta))) val = 2;
+ else val = 1 + (long long)(n * pow(eta*u - eta + 1., alpha));
+
+ return(val);
+}
+
+/*
+long double zipf(double theta, long n) {
+
+ // this draws values from the zipf distribution
+ // this is mainly useful for test generation purposes
+ // n is range, theta is skewness parameter
+ // theta = 0 gives uniform dbn,
+ // theta > 1 gives highly skewed dbn.
+
+ double alpha;
+ double zetan;
+ double eta;
+ double u;
+ double uz;
+ long double val;
+
+ // randinit must be called before entering this procedure for
+ // the first time since it uses the random generators
+
+ alpha = 1. / (1. - theta);
+ zetan = zeta(n, theta);
+ eta = (1. - pow(2./n, 1. - theta)) / (1. - zeta(2.,theta)/zetan);
+
+ u = randomfl();
+ uz = u * zetan;
+ if (uz < 1.) val = 1;
+ else if (uz < (1. + pow(0.5, theta))) val = 2;
+ else val = 1 + (long long)(n * pow(eta*u - eta + 1., alpha));
+
+ return(val);
+}
+*/
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/debian/prng.h ntop-4.99.3+ndpi5517+dfsg2/debian/prng.h
--- ntop-4.99.3+ndpi5517+dfsg1/debian/prng.h 1969-12-31 16:00:00.000000000 -0800
+++ ntop-4.99.3+ndpi5517+dfsg2/debian/prng.h 2012-11-30 00:37:15.000000000 -0800
@@ -0,0 +1,53 @@
+// Probabilistic Random Number Generators
+// Collected from various sources by Graham Cormode, 2000-2003
+//
+
+#include <math.h>
+
+#ifndef _PRNG
+
+#define MOD 2147483647
+#define HL 31
+
+extern long hash31(long long, long long, long long);
+extern long fourwise(long long, long long, long long, long long, long long);
+
+#define KK 17
+#define NTAB 32
+
+typedef struct prng_type{
+ int usenric; // which prng to use
+ float scale; /* 2^(- integer size) */
+ long floatidum;
+ long intidum; // needed to keep track of where we are in the
+ // nric random number generators
+ long iy;
+ long iv[NTAB];
+ /* global variables */
+ unsigned long randbuffer[KK]; /* history buffer */
+ int r_p1, r_p2; /* indexes into history buffer */
+ int iset;
+ double gset;
+} prng_type;
+
+#define _PRNG 1
+
+#endif
+
+extern long prng_int(prng_type *);
+extern float prng_float(prng_type *);
+extern prng_type * prng_Init(long, int);
+extern void prng_Destroy(prng_type * prng);
+void prng_Reseed(prng_type *, long);
+
+//extern long double zipf(double, long) ;
+extern double fastzipf(double, long, double, prng_type *);
+extern double zeta(long, double);
+extern double prng_normal(prng_type * prng);
+extern double prng_stable(prng_type * prng, double);
+
+//extern double stable(double); // stable distributions
+//extern double stabledevd(float) ;
+//extern long double stabledevl(float) ;
+//extern double altstab(double);
+
diff -Nru ntop-4.99.3+ndpi5517+dfsg1/prng.c ntop-4.99.3+ndpi5517+dfsg2/prng.c
--- ntop-4.99.3+ndpi5517+dfsg1/prng.c 2012-02-21 15:51:02.000000000 -0800
+++ ntop-4.99.3+ndpi5517+dfsg2/prng.c 2012-11-30 00:16:18.000000000 -0800
@@ -1,3 +1,24 @@
+/*********************************************************************
+
+prng: pseudo-random number generators and hashing routines
+
+G. Cormode 2003-2012
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License along
+ with this program; if not, write to the Free Software Foundation, Inc.,
+ 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+************************************************************************/
+
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
@@ -9,7 +30,7 @@
{
long long result;
- long lresult;
+ long lresult;
// return a hash of x using a and b mod (2^31 - 1)
// may need to do another mod afterwards, or drop high bits
@@ -19,8 +40,8 @@
// result = ((long long) a)*((long long) x)+((long long) b);
result=(a * x) + b;
result = ((result >> HL) + result) & MOD;
- lresult=(long) result;
-
+ lresult=(long) result;
+
return(lresult);
}
@@ -28,11 +49,11 @@
{
long long result;
long lresult;
-
+
// returns values that are 4-wise independent by repeated calls
- // to the pairwise indpendent routine.
+ // to the pairwise indpendent routine.
- result = hash31(hash31(hash31(x,a,b),x,c),x,d);
+ result = hash31(hash31(hash31(a,b,x),c,x),d,x);
lresult = (long) result;
return lresult;
}
@@ -45,7 +66,7 @@
// There are *THREE* alternate implementations of PRNGs here.
// One taken from Numerical Recipes in C, the second from www.agner.org
// The third is an internal C random library, srand
-// The variable usenric controls which one is used: pick one
+// The variable usenric controls which one is used: pick one
// and stick with it, switching between the two will give unpredictable
// results. This is controlled by the randinit procedure, call it with
// usenric == 1 to use the Numerical Recipes gens
@@ -66,7 +87,7 @@
#define RNMX (1.0-EPS)
float ran1(prng_type * prng) {
-
+
// A Random Number Generator that picks a uniform [0,1] random number
// From Numerical Recipes, page 280
// Should be called with a NEGATIVE value of idum to initialize
@@ -75,7 +96,7 @@
int j;
long k;
float temp;
-
+
if (prng->floatidum <= 0 || !prng->iy) {
if (-(prng->floatidum) < 1) prng->floatidum=1;
else prng->floatidum = -(prng->floatidum);
@@ -96,9 +117,9 @@
if ((temp=AM*prng->iy) > RNMX) return RNMX;
else return temp;
}
-
+
long ran2(prng_type * prng) {
-
+
// A Random Number Generator that picks a uniform random number
// from the range of long integers.
// From Numerical Recipes, page 280
@@ -109,7 +130,7 @@
int j;
long k;
-
+
if (prng->intidum <= 0 || !prng->iy) {
if (-(prng->intidum) < 1) prng->intidum=1;
else prng->intidum = -(prng->intidum);
@@ -174,7 +195,7 @@
/* generate next random number */
- x = prng->randbuffer[prng->r_p1] = rotl(prng->randbuffer[prng->r_p2], R1)
+ x = prng->randbuffer[prng->r_p1] = rotl(prng->randbuffer[prng->r_p2], R1)
+ rotl(prng->randbuffer[prng->r_p1], R2);
/* rotate list pointers */
if (--prng->r_p1 < 0) prng->r_p1 = KK - 1;
@@ -225,7 +246,7 @@
{
case 1 : response=(ran2(prng)); break;
case 2 : response=(ran3(prng)); break;
- // case 3 : response=(lrand48()); break;
+ case 3 : response=(lrand48()); break;
}
return response;
}
@@ -233,7 +254,7 @@
float prng_float(prng_type * prng) {
- // returns a pseudo-random float in the range [0.0,1.0].
+ // returns a pseudo-random float in the range [0.0,1.0].
// Initialise the generator before use!
float result=0;
@@ -241,7 +262,7 @@
{
case 1 : result=(ran1(prng)); break;
case 2 : result=(ran4(prng)); break;
- // case 3 : result=(drand48()); break;
+ case 3 : result=(drand48()); break;
}
return result;
}
@@ -249,25 +270,25 @@
prng_type * prng_Init(long seed, int nric) {
// Initialise the random number generators. nric determines
- // which algorithm to use, 1 for Numerical Recipes in C,
- // 0 for the other one.
+ // which algorithm to use, 1 for Numerical Recipes in C,
+ // 0 for the other one.
prng_type * result;
result=(prng_type *) calloc(1,sizeof(prng_type));
result->iy=0;
- result->usenric=nric;
+ result->usenric=nric;
result->floatidum=-1;
result->intidum=-1;
result->iset=0;
// set a global variable to record which algorithm to use
switch (nric)
- {
- case 2 :
+ {
+ case 2 :
RanrotAInit(result,seed);
break;
- case 1 :
- if (seed>0) {
+ case 1 :
+ if (seed>0) {
// to initialise the NRiC PRNGs, call it with a negative value
// so make sure it gets a negative value!
result->floatidum = -(seed); result->intidum = -(seed);
@@ -275,26 +296,26 @@
result->floatidum=seed; result->intidum=seed;
}
break;
- case 3 :
- srand(seed);
+ case 3 :
+ srand48(seed);
break;
}
prng_float(result);
prng_int(result);
// call the routines to actually initialise them
- return(result);
+ return(result);
}
void prng_Reseed(prng_type * prng, long seed)
{
switch (prng->usenric)
- {
- case 2 :
+ {
+ case 2 :
RanrotAInit(prng,seed);
break;
- case 1 :
- if (seed>0) {
+ case 1 :
+ if (seed>0) {
// to initialise the NRiC PRNGs, call it with a negative value
// so make sure it gets a negative value!
prng->floatidum = -(seed); prng->intidum = -(seed);
@@ -302,16 +323,15 @@
prng->floatidum=seed; prng->intidum=seed;
}
break;
- case 3 :
- srand(seed);
+ case 3 :
+ srand48(seed);
break;
}
-}
+}
void prng_Destroy(prng_type * prng)
-{
+{
free(prng);
- prng=NULL;
}
/**********************************************************************/
@@ -333,7 +353,7 @@
v2 = 2.0*prng_float(prng)-1.0;
rsq=v1*v1+v2*v2;
} while (rsq >= 1.0 || rsq == 0.0);
-
+
fac = sqrt((double) -2.0*log((double)rsq)/rsq);
prng->gset=v1*fac;
prng->iset=1;
@@ -349,7 +369,7 @@
// From 'stable distributions', John Nolan, manuscript, p24
// we set beta = 0 by analogy with the normal and cauchy case
- // identical to the above routine, but returns a double instead
+ // identical to the above routine, but returns a double instead
// of a long double (you'll see this a lot...)
double theta, W, holder, left, right;
@@ -380,10 +400,10 @@
}
-double prng_altstab(prng_type * prng, double p)
+double prng_altstab(prng_type * prng, double p)
{
double u,v,result;
-
+
u=prng_float(prng);
v=prng_float(prng);
result=pow(u,p);
@@ -402,7 +422,7 @@
z=gasdev();
return (1.0/(z*z));
- }
+ }
*/
@@ -416,7 +436,7 @@
// the first time since it uses the random generators
- if (alpha==2.0)
+ if (alpha==2.0)
return(prng_normal(prng));
else if (alpha==1.0)
return(prng_cauchy(prng));
@@ -425,7 +445,7 @@
else return (prng_stabledbn(prng,alpha));
}
-double zeta(long n, double theta)
+double zeta(long n, double theta)
{
// the zeta function, used by the below zipf function
@@ -434,7 +454,7 @@
int i;
double ans=0.0;
-
+
for (i=1; i <= n; i++)
ans += pow(1./(double)i, theta);
return(ans);
@@ -444,9 +464,9 @@
// this draws values from the zipf distribution
// this is mainly useful for test generation purposes
- // n is range, theta is skewness parameter
+ // n is range, theta is skewness parameter
// theta = 0 gives uniform dbn,
- // theta > 1 gives highly skewed dbn.
+ // theta > 1 gives highly skewed dbn.
// original code due to Flip Korn, used with permission
double alpha;
@@ -459,7 +479,8 @@
// the first time since it uses the random generators
alpha = 1. / (1. - theta);
- eta = (1. - pow(2./n, 1. - theta)) / (1. - zeta(2.,theta)/zetan);
+ eta = (1. - pow(2./((double) n), 1. - theta))
+ / (1. - zeta(2,theta)/zetan);
u = prng_float(prng);
uz = u * zetan;
@@ -475,9 +496,9 @@
// this draws values from the zipf distribution
// this is mainly useful for test generation purposes
- // n is range, theta is skewness parameter
+ // n is range, theta is skewness parameter
// theta = 0 gives uniform dbn,
- // theta > 1 gives highly skewed dbn.
+ // theta > 1 gives highly skewed dbn.
double alpha;
double zetan;
--- End Message ---