00001 #include <vector>
00002 #include <algorithm>
00003 #include "CommonTools/Statistics/interface/SequentialPartitionGenerator.h"
00004
00005 SequentialPartitionGenerator::SequentialPartitionGenerator( int n, int k, int pmin ) :
00006 the_n(n), the_k(k), the_pmin(pmin), the_pmax(n) {}
00007
00008 SequentialPartitionGenerator::SequentialPartitionGenerator(int n, int k, int pmin, int pmax) :
00009 the_n(n), the_k(k), the_pmin(pmin), the_pmax(pmax) {}
00010
00011 SequentialPartitionGenerator::Partition SequentialPartitionGenerator::next_partition()
00012 {
00013 bool done=next_part(the_part);
00014 if (done)
00015 {
00016 return the_part;
00017 };
00018 SequentialPartitionGenerator::Partition empty;
00019 return empty;
00020 }
00021
00022 bool SequentialPartitionGenerator::first_part(
00023 SequentialPartitionGenerator::Partition & p, int k, int n, int pmin, int pmax ) const
00024 {
00025 n_first++;
00026 bool done=false;
00027 switch (k) {
00028 case 1:
00029 p[0]=std::min(pmax,n);
00030 return (n<=pmax && p[0]>=pmin);
00031 case 2:
00032 for (int i=std::min(pmax,n-1);i>=pmin;i--) {
00033 if (done=(i>=n-i && n-i>=pmin)) {p[0]=i;p[1]=n-i;}
00034 return done;
00035 }
00036 default:
00037 Partition pp(p.begin()+1,p.end());
00038 for (int i=std::min(pmax,n-k+1);i>=pmin;i--) {
00039 p[0]=i;
00040 done=this->first_part(pp,k-1,n-p[0],pmin,p[0]);
00041 if (done) {copy(pp.begin(),pp.end(),p.begin()+1);}
00042 return done;
00043 }
00044 }
00045 return done;
00046 }
00047
00048 bool SequentialPartitionGenerator::next_part(
00049 SequentialPartitionGenerator::Partition & p ) const
00050 {
00051 n_next++;
00052 bool done=false;
00053 int k=p.size();
00054 switch (k) {
00055 case 0:
00056 p.insert(p.begin(),the_k,0);
00057 return this->first_part(p,the_k,the_n,the_pmin,the_pmax);
00058 case 1:
00059 return false;
00060 default:
00061 int n=0;
00062 for (int i=0;i<k;i++) n=n+p[i];
00063 SequentialPartitionGenerator::Partition pp(p.begin()+1,p.end());
00064 done = (pp.size()>1 ? this->next_part(pp) : false);
00065 if (done)
00066 {
00067 copy(pp.begin(),pp.end(),p.begin()+1);
00068 } else {
00069 done = (p[0]==1 ? false :
00070 this->first_part(pp,k-1,n-p[0]+1,the_pmin,p[0]-1));
00071 if (done) { --p[0];copy(pp.begin(),pp.end(),p.begin()+1); }
00072 };
00073 return done;
00074 };
00075 return done;
00076 }