CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_0/src/CommonTools/Statistics/interface/SequentialCombinationGenerator.h

Go to the documentation of this file.
00001 #ifndef SequentialCombinationGenerator_H
00002 #define SequentialCombinationGenerator_H
00003 
00004 #include "CommonTools/Statistics/interface/SequentialPartitionGenerator.h"
00005 #include <vector>
00006 #include <list>
00007 #include <algorithm>
00008 #include <numeric>
00009 #include <functional>
00010 #include <iostream>
00011 
00021 template<class T> class SequentialCombinationGenerator {
00022 public:
00023   typedef SequentialPartitionGenerator::Partition Partition;
00024   typedef std::vector<T> Collection;
00025   typedef std::vector<Collection> Combination;
00026   typedef std::vector<int> Vecint;
00027 
00028   SequentialCombinationGenerator ( Partition & part ) :
00029       the_part ( part ), the_k ( part.size() ),
00030       the_n ( accumulate ( part.begin(), part.end(), 0 ) )
00031   {
00032     sort(the_part.begin(),the_part.end());
00033     the_comb.reserve(the_n);
00034   };
00035 
00039   Combination next_combination ( Collection& coll )
00040   {
00041     dbg=false;
00042     Combination comb;
00043     Vecint newcomb=next_combi(the_comb,the_n,the_part);
00044     the_comb=newcomb;
00045     if (newcomb.size()==0) {
00046       Combination empty;
00047       return empty;
00048     }
00049     int i=0;
00050     for (int j=0;j<the_k;j++) 
00051     {
00052       Collection temp;
00053       for (int l=0;l<the_part[j];l++) 
00054       {
00055         temp.push_back(coll[the_comb[i]-1]);
00056         i++;
00057       }
00058       comb.push_back(temp);
00059     }
00060     return comb;
00061   };
00062 
00063 private:
00064   Vecint next_combi( Vecint & cold, int n, const Partition & p )
00065   {
00066     Vecint empty;
00067     if (cold.size()==0) { // first entry, initialize
00068       cold.reserve(n);
00069       for (int i=0;i<n;cold.push_back(++i));
00070       return cold;
00071     }
00072     int k=p.size();
00073     if (k==1) return empty;
00074     Vecint cnew(cold);
00075     int n1=n-p[0];
00076     Vecint cold1(cold.begin()+p[0],cold.end());
00077     Vecint cnew1(n1);
00078     Partition p1(p.begin()+1,p.end());
00079     cnew1=next_combi(cold1,n1,p1);
00080     if (cnew1.size()!=0) {
00081       copy(cnew1.begin(),cnew1.end(),cnew.begin()+p[0]);
00082       return cnew;
00083     }
00084     Vecint cold2(cold.begin(),cold.begin()+p[0]);
00085     Vecint cnew2(p[0]);
00086     sort(cold.begin(),cold.end());
00087     sort(cold2.begin(),cold2.end());
00088     cnew2=next_subset(cold,cold2);
00089     if (cnew2.size()==0) return empty;
00090     copy(cnew2.begin(),cnew2.begin()+p[0],cnew.begin());
00091     Vecint ss(n1);
00092     set_difference(cold.begin(),cold.end(),
00093                    cnew2.begin(),cnew2.end(),ss.begin());
00094     int ip=p[0];
00095     for (int i=1;i<k;i++) {
00096       if (p[i]!=p[i-1]) {
00097         copy(ss.begin(),ss.end(),&cnew[ip]);
00098         return cnew;
00099       }
00100       int mincnew2=cnew2[0];
00101       if (ss[n1-1]<mincnew2) return empty;
00102       Vecint::iterator j1=find_if(ss.begin(),ss.end(),bind2nd(std::greater<int>(),mincnew2));
00103       if (ss.end()-j1 < p[i]) return empty;
00104       Vecint sss(j1,ss.end());
00105       for (int j=0;j<p[i];cnew[ip+j]=cnew2[j]=sss[j++]);
00106       int n2=ss.size()-cnew2.size();
00107       if (n2==0) return cnew;
00108       Vecint s(n2);
00109       set_difference(ss.begin(),ss.end(),cnew2.begin(),cnew2.end(),s.begin());
00110       ss=s;
00111       ip+=p[i];
00112     }
00113    
00114    return empty; 
00115  
00116   };
00117 
00118   Vecint next_subset(Vecint g, Vecint c)
00119   {
00120     Vecint empty;
00121     int n=g.size();
00122     int k=c.size();
00123     typename Vecint::iterator ind;
00124     for (int i=k-1;i>=0;i--) {
00125       if (c[i]<g[n-k+i]) {
00126         
00127         
00128 //      ind=find(&g[0],&g[n-k+i],c[i])+1;
00129         
00130         Vecint::iterator g2 = g.begin();
00131         advance(g2,n-k+i);  
00132         
00133         ind=find(g.begin(),g2,c[i]);
00134         
00135         ind++;
00136         
00137         
00138         
00139         copy(ind,ind+k-i,&c[i]);
00140         return c;
00141       }
00142     }
00143     return empty;
00144   };
00145 
00146   void vecprint(const Vecint & v) const
00147   {
00148     int n=v.size();
00149     for (int i=0;i<n;std::cout << v[i++]);
00150     std::cout << std::endl;
00151   };
00152 
00153 private:
00154   int the_n;
00155   int the_k;
00156   Partition the_part;
00157   mutable Vecint the_comb;
00158   mutable bool dbg;
00159 };
00160 
00161 #endif