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