Index: /trunk/Cosmo/RadioBeam/anapkn.pic
===================================================================
--- /trunk/Cosmo/RadioBeam/anapkn.pic	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/anapkn.pic	(revision 3947)
@@ -25,5 +25,5 @@
 z = 0.5
 # On considere une fraction de HI (HI/Baryon) = 0.02
-fHI = 0.020/0.1 
+gHI = 0.02/0.01
 
 ####################################################################
@@ -31,10 +31,10 @@
 defscript convpk2t21 '  --> calcul le coefficient de conversion P(k) en P(k)_temperature mK^2'
   cgeom = (1+$z)*(1+$z)/sqrt($OL+pow((1+$z),3.)*$Om)
-  cct21 = 0.57*$cgeom*$fHI
+  cct21 = 0.054*$cgeom*$gHI
   cct21 = $cct21*$cct21
   echo "----  cct21 (NoFactCroiss)= $cct21"
 #  On approxime le facteur de croissance lineaire au carre par 
 #  [(1+0.5)/(1+z)]^2  --> erreur ~ 10%
-  cct21 = $cct21*pow(1.5/(1+$z),2)
+#  cct21 = $cct21*pow(1.5/(1+$z),2)
   echo "---- Conv P(k)->Temperature mK^2 : cct21= $cct21"
 endscript
Index: /trunk/Cosmo/RadioBeam/interfconfigs.cc
===================================================================
--- /trunk/Cosmo/RadioBeam/interfconfigs.cc	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/interfconfigs.cc	(revision 3947)
@@ -168,6 +168,6 @@
 
  return vd;
-
-}
+}
+
 /* --Fonction -- */
 vector<Dish> CreateConfigD(double Ddish, double Eta)
@@ -243,6 +243,105 @@
 
  cout << ">>>CreateConfigD(" << Ddish << "," << Eta << ") ---> NDishes=" << vd.size() << endl;
-
  return vd;
+}
+
+/* --Fonction -- */
+vector<Dish> CreateConfigNancay12(double Ddish, double Eta)
+{
+  EnumeratedSequence esx,esy;
+  esx = 0,1,3;
+  esy = 0,1,2,3;
+
+  vector<Dish> vd;
+  int cnt=0;
+
+  for(int ix=0; ix<esx.Size(); ix++) {
+    for(int iy=0; iy<esy.Size(); iy++) {
+      cnt++; 
+      vd.push_back(Dish(cnt, ((double)esx.Value(ix))*Ddish,((double)esy.Value(iy))*Ddish,Eta*Ddish));
+    }
+  }
+  cout << ">>>CreateConfigNancay12(" << Ddish << "," << Eta << ") ---> NDishes=" << vd.size() << endl;
+  return vd;
+}
+
+/* --Fonction -- */
+vector<Dish> CreateConfigNancay24(double Ddish, double Eta)
+{
+  EnumeratedSequence esx,esy;
+  esx = 0,1,3,5;
+  esy = 0,1,2,3,4,5;
+
+  vector<Dish> vd;
+  int cnt=0;
+
+  for(int ix=0; ix<esx.Size(); ix++) {
+    for(int iy=0; iy<esy.Size(); iy++) {
+      cnt++; 
+      vd.push_back(Dish(cnt, ((double)esx.Value(ix))*Ddish,((double)esy.Value(iy))*Ddish,Eta*Ddish));
+    }
+  }
+  cout << ">>>CreateConfigNancay24(" << Ddish << "," << Eta << ") ---> NDishes=" << vd.size() << endl;
+  return vd;
+}
+
+/* --Fonction -- */
+vector<Dish> CreateConfigNancay36(double Ddish, double Eta)
+{
+  EnumeratedSequence esx,esy;
+  esx = 0,1,2,3,4,5; 
+  esy = 0,1,2,3,4,5; 
+
+  vector<Dish> vd;
+  int cnt=0;
+
+  for(int ix=0; ix<esx.Size(); ix++) {
+    for(int iy=0; iy<esy.Size(); iy++) {
+      cnt++; 
+      vd.push_back(Dish(cnt, ((double)esx.Value(ix))*Ddish,((double)esy.Value(iy))*Ddish,Eta*Ddish));
+    }
+  }
+  cout << ">>>CreateConfigNancay36(" << Ddish << "," << Eta << ") ---> NDishes=" << vd.size() << endl;
+  return vd;
+}
+
+/* --Fonction -- */
+vector<Dish> CreateConfigNancay40(double Ddish, double Eta)
+{
+  EnumeratedSequence esx,esy;
+  esx = 0,1,3,5,7;
+  esy = 0,1,2,3,4,5,6,7;
+
+  vector<Dish> vd;
+  int cnt=0;
+
+  for(int ix=0; ix<esx.Size(); ix++) {
+    for(int iy=0; iy<esy.Size(); iy++) {
+      cnt++; 
+      vd.push_back(Dish(cnt, ((double)esx.Value(ix))*Ddish,((double)esy.Value(iy))*Ddish,Eta*Ddish));
+    }
+  }
+  cout << ">>>CreateConfigNancay40(" << Ddish << "," << Eta << ") ---> NDishes=" << vd.size() << endl;
+  return vd;
+}
+
+/* --Fonction -- */
+vector<Dish> CreateConfigNancay128(double Ddish, double Eta)
+{
+  EnumeratedSequence esx,esy;
+  esx = 0,1,3,4,7,10,13,15; 
+  esy = 0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15;
+
+  vector<Dish> vd;
+  int cnt=0;
+
+  for(int ix=0; ix<esx.Size(); ix++) {
+    for(int iy=0; iy<esy.Size(); iy++) {
+      cnt++; 
+      vd.push_back(Dish(cnt, ((double)esx.Value(ix))*Ddish,((double)esy.Value(iy))*Ddish,Eta*Ddish));
+    }
+  }
+  cout << ">>>CreateConfigNancay128(" << Ddish << "," << Eta << ") ---> NDishes=" << vd.size() << endl;
+  return vd;
 }
 
Index: /trunk/Cosmo/RadioBeam/interfconfigs.h
===================================================================
--- /trunk/Cosmo/RadioBeam/interfconfigs.h	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/interfconfigs.h	(revision 3947)
@@ -26,4 +26,15 @@
 vector<Dish> CreateConfigD(double Ddish=5., double Eta=0.9);
 
+// Configuration a 12 dishs - arrange en cylindres (4+4+4) 
+vector<Dish> CreateConfigNancay12(double Ddish=5., double Eta=0.9);
+// Configuration a 24 dishs - arrange en cylindres (6+6+0+6+0+6) 
+vector<Dish> CreateConfigNancay24(double Ddish=5., double Eta=0.9);
+// Configuration a 36 cylindres (6+6+6+6+6+6) 
+vector<Dish> CreateConfigNancay36(double Ddish=5., double Eta=0.9);
+// Configuration a 40 cylindres (8+8+0+8+0+8+0+8) 
+vector<Dish> CreateConfigNancay40(double Ddish=5., double Eta=0.9);
+
+// Configuration a 128 dishs - arrange en 8 cylindres de 15 dishs (CC0CC00C00C00C0C)
+vector<Dish> CreateConfigNancay128(double Ddish=5., double Eta=0.9);
 
 // Filled cylinder configuration 
Index: /trunk/Cosmo/RadioBeam/mdish.cc
===================================================================
--- /trunk/Cosmo/RadioBeam/mdish.cc	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/mdish.cc	(revision 3947)
@@ -65,4 +65,48 @@
 }
 
+HProf Four2DResponse::GetProjNoiseLevel(int nbin, bool fgnorm1)
+{
+  Histo2D h2w = GetResponse(2*nbin, 2*nbin);
+  r_8 vmin=h2w.VMin();
+  r_8 vmax=h2w.VMax();
+  double seuil=vmax/10000.;
+  if (seuil<1.e-6) seuil=1.e-6;
+  r_8 facnorm=((fgnorm1)?vmax:1.);
+  cout << " Four2DResponse::GetProjNoiseLevel Min,Max=" << vmin << " , " << vmax 
+       << " facnorm=" << facnorm << " seuil=" << seuil << endl;  
+  double kmax=2.*M_PI*D();
+  HProf hp(0,kmax,nbin);
+  double x,y;
+  for(sa_size_t j=0; j<h2w.NBinY(); j++) {
+    for(sa_size_t i=0; i<h2w.NBinX(); i++) {
+      h2w.BinCenter(i,j,x,y);
+      double yw=h2w(i,j);
+      if (yw<seuil) continue;
+      hp.Add(sqrt(x*x+y*y),facnorm/yw);
+    }
+  }
+  return hp;
+}
+
+HProf Four2DResponse::GetProjResponse(int nbin, bool fgnorm1)
+{
+  Histo2D h2w = GetResponse(2*nbin, 2*nbin);
+  r_8 vmin=h2w.VMin();
+  r_8 vmax=h2w.VMax();
+  r_8 facnorm=((fgnorm1)?(1./vmax):1.);
+  cout << " Four2DResponse::GetProjResponse Min,Max=" << vmin << " , " << vmax 
+       << " facnorm=" << facnorm << endl; 
+  double kmax=2.*M_PI*D();
+  HProf hp(0,kmax,nbin);
+  double x,y;
+  for(sa_size_t j=0; j<h2w.NBinY(); j++) {
+    for(sa_size_t i=0; i<h2w.NBinX(); i++) {
+      h2w.BinCenter(i,j,x,y);
+      hp.Add(sqrt(x*x+y*y),h2w(i,j)*facnorm);
+    }
+  }
+  return hp;
+}
+
 //---------------------------------------------------------------
 // -- Four2DRespTable : Reponse tabulee instrumentale ds le plan k_x,k_y (angles theta,phi) 
@@ -185,10 +229,11 @@
   SetBeamNSamples();
   SetPrtLevel();
+  fgcomputedone_=false;
   mcnt_=0;
 }
 
-Histo2D MultiDish::GetResponse()
-{
-  cout << " MultiDish::GetResponse() - NDishes=" << dishes_.size() << " nx=" << nx_ << " ny=" << ny_ << endl;
+void MultiDish::ComputeResponse()
+{
+  cout << " MultiDish::ComputeResponse() - NDishes=" << dishes_.size() << " nx=" << nx_ << " ny=" << ny_ << endl;
   double kmx = 1.2*DeuxPI*dmax_/lambda_;
   double dkmx = kmx/(double)nx_;
@@ -198,4 +243,5 @@
   h2w_.Define(-kmxx,kmxx,2*nx_+1,-kmxy,kmxy,2*ny_+1);
   h2w_.SetZeroBin(0.,0.);
+  kmax_=kmx;
 
   double dold = dishes_[0].Diameter()/lambda_;
@@ -210,5 +256,5 @@
   double dtet = thetamax_/(double)ntet_;
   double dphi = phimax_/(double)nphi_;
-  cout << " MultiDish::GetResponse() - ThetaMax=" << thetamax_ << " NTheta=" << ntet_ 
+  cout << " MultiDish::ComputeResponse() - ThetaMax=" << thetamax_ << " NTheta=" << ntet_ 
        << " PhiMax=" <<  phimax_ << " NPhi=" << nphi_ << endl;
 
@@ -224,6 +270,20 @@
       sumw += CumulResp(rd, theta, -(phi+M_PI));
     }
-    if (prtlev_>0) cout << "  MultiDish::GetResponse() done ktheta=" << kt << " / MaxNTheta= " << ntet_ << endl; 
-  }
+    if (prtlev_>0) 
+      cout << "  MultiDish::ComputeResponse() done ktheta=" << kt << " / MaxNTheta= " 
+	   << ntet_ << endl; 
+  }
+  r_8 rmin,rmax;
+  h2w_.GetMinMax(rmin,rmax);
+  cout << "  MultiDish::ComputeResponse() finished : Rep_min,max=" << rmin << "," << rmax << " sumW0=" 
+       << sumw << " ?=" << h2w_.SumWBinZero() << endl;
+  fgcomputedone_=true;
+  return;
+}
+
+Histo2D MultiDish::GetResponse()
+{
+  if (!fgcomputedone_) ComputeResponse();
+
   double kx1 = DeuxPI*(dishes_[0].DiameterX())/lambda_;
   double ky1 = DeuxPI*(dishes_[0].DiameterY())/lambda_;
@@ -236,4 +296,5 @@
     ib=jb=0;
   }
+  double sumw=h2w_.SumWBinZero();
   double vmax=h2.VMax();
   cout << " MultiDish::GetResponse[1] VMin=" << h2.VMin() << " VMax= " << vmax  
@@ -257,4 +318,45 @@
   return h2;
 }
+
+HProf MultiDish::GetProjNoiseLevel(int nbin, bool fgnorm1)
+{
+  r_8 vmin,vmax;
+  h2w_.GetMinMax(vmin,vmax);
+  double seuil=vmax/10000.;
+  if (seuil<1.e-6) seuil=1.e-6;
+  r_8 facnorm=((fgnorm1)?vmax:1.);
+  cout << " MultiDish::GetProjNoiseLevel Min,Max=" << vmin << " , " << vmax 
+       << " facnorm=" << facnorm << " seuil=" << seuil << endl;  
+  HProf hp(0,kmax_,nbin);
+  for(sa_size_t j=0; j<h2w_.NBinY(); j++) {
+    double y=h2w_.Y(j);
+    for(sa_size_t i=0; i<h2w_.NBinX(); i++) {
+      double x=h2w_.X(i);
+      double yw=h2w_(i,j);
+      if (yw<seuil) continue;
+      hp.Add(sqrt(x*x+y*y),facnorm/yw);
+    }
+  }
+  return hp;
+}
+
+HProf MultiDish::GetProjResponse(int nbin, bool fgnorm1)
+{
+  r_8 vmin,vmax;
+  h2w_.GetMinMax(vmin,vmax);
+  r_8 facnorm=((fgnorm1)?(1./vmax):1.);
+  cout << " MultiDish::GetProjResponse Min,Max=" << vmin << " , " << vmax 
+       << " facnorm=" << facnorm << endl;  
+  HProf hp(0,kmax_,nbin);
+  for(sa_size_t j=0; j<h2w_.NBinY(); j++) {
+    double y=h2w_.Y(j);
+    for(sa_size_t i=0; i<h2w_.NBinX(); i++) {
+      double x=h2w_.X(i);
+      hp.Add(sqrt(x*x+y*y),h2w_(i,j)*facnorm);
+    }
+  }
+  return hp;
+}
+
 
 Histo2D MultiDish::PosDist(int nx, int ny, double dmax)
@@ -312,4 +414,5 @@
       double kx0 = DeuxPI*(dishes_[i].X-dishes_[j].X)/lambda_;
       double ky0 = DeuxPI*(dishes_[i].Y-dishes_[j].Y)/lambda_;
+      double pgain=dishes_[i].Gain()*dishes_[j].Gain();
       rot.Do(kx0, ky0);
       //    if (kx0<0) kx0=-kx0;
@@ -321,14 +424,14 @@
 	  double y = jy*dy;
 	  if ((ix>0)&&(jy>0)) {
-	    sumw += AddToHisto(kx0, ky0, x, y, rd(x,y), fgfh);
-	    if (j!=i) sumw += AddToHisto(-kx0, -ky0, x, y, rd(x,y), fgfh);
+	    sumw += AddToHisto(kx0, ky0, x, y, rd(x,y)*pgain, fgfh);
+	    if (j!=i) sumw += AddToHisto(-kx0, -ky0, x, y, rd(x,y)*pgain, fgfh);
 	  }
 	  else {
 	    if ((ix==0)&&(jy==0)) {
-	      sumw += h2w_.Add(kx0, ky0, rd(0.,0.), fgfh);
-	      if (j!=i) sumw += h2w_.Add(-kx0, -ky0, rd(0.,0.), fgfh);
+	      sumw += h2w_.Add(kx0, ky0, rd(0.,0.)*pgain, fgfh);
+	      if (j!=i) sumw += h2w_.Add(-kx0, -ky0, rd(0.,0.)*pgain, fgfh);
 	    }
 	    else {
-	      double w = rd(x,y);
+	      double w = rd(x,y)*pgain;
 	      if (ix==0) {
 		sumw += h2w_.Add(kx0, ky0+y, w, fgfh);
Index: /trunk/Cosmo/RadioBeam/mdish.h
===================================================================
--- /trunk/Cosmo/RadioBeam/mdish.h	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/mdish.h	(revision 3947)
@@ -42,4 +42,7 @@
   inline void setLambda(double lambda=1.)
   { lambda_ = lambda;   lambda_ratio_ = lambda_/lambdaref_; }   
+
+  inline double getLambdaRef()   { return lambdaref_; }
+  inline double getLambda()  { return lambda_; } 
   
   // Return the 2D response for wave vector (kx,ky)
@@ -48,4 +51,9 @@
   { return Value(kx, ky); } 
   virtual Histo2D GetResponse(int nx=255, int ny=255);  
+  // Retourne le niveau moyen du bruit projete 1D en fonction (sqrt(u^2+v^2)
+  HProf GetProjNoiseLevel(int nbin=128, bool fgnorm1=true);
+  // Retourne la reponse moyenne projetee 1D en fonction (sqrt(u^2+v^2)
+  HProf GetProjResponse(int nbin=128, bool fgnorm1=true);
+
   inline double D() { return dx_; } ; 
   inline double Dx() { return dx_; } ; 
@@ -95,10 +103,11 @@
   // Circular dish 
   Dish(int id, double x, double y, double diam)
-    :  id_(id), X(x), Y(y), D(diam), Dx(D), Dy(D), fgcirc_(true)   {   } 
+    :  id_(id), X(x), Y(y), D(diam), Dx(D), Dy(D), fgcirc_(true), gain_(1.)   {   } 
   // Receiver with rectangular type answer in kx,ky plane 
   Dish(int id, double x, double y, double dx, double dy)
-    :  id_(id), X(x), Y(y), D(sqrt(dx*dy)), Dx(dx), Dy(dy), fgcirc_(false)   {   } 
+    :  id_(id), X(x), Y(y), D(sqrt(dx*dy)), Dx(dx), Dy(dy), fgcirc_(false), gain_(1.)    {   } 
   Dish(Dish const& a) 
-    :  id_(a.id_), X(a.X), Y(a.Y), D(a.D), Dx(a.Dx), Dy(a.Dy), fgcirc_(a.fgcirc_)     {   } 
+    :  id_(a.id_), X(a.X), Y(a.Y), D(a.D), Dx(a.Dx), Dy(a.Dy), fgcirc_(a.fgcirc_), gain_(a.gain_)  {   } 
+  inline void setGain(double gain) { gain_=gain; return; }
   inline bool isCircular() { return fgcirc_; }
   inline int ReflectorId() { return id_; }
@@ -106,4 +115,5 @@
   inline double DiameterX() { return Dx; }
   inline double DiameterY() { return Dy; }
+  inline double Gain() { return gain_; }
 
   int id_;   // numero de reflecteur
@@ -111,4 +121,5 @@
   double Dx, Dy;
   bool fgcirc_;  // false -> rectangular dish 
+  double gain_;
 };
 
@@ -132,5 +143,13 @@
   { beamnx_=nx; beamny_=ny; }
 
+  // Calcul la reponse ds le plan 2D  (u,v) = (kx,ky)
+  void ComputeResponse();
+  // Retourne la reponse 2D ds le plan (u,v) = (kx,ky) sous forme d'histo 2D 
   Histo2D GetResponse();
+
+  // Retourne le niveau moyen du bruit projete 1D en fonction (sqrt(u^2+v^2)
+  HProf GetProjNoiseLevel(int nbin=128, bool fgnorm1=true);
+  // Retourne la reponse moyenne projetee 1D en fonction (sqrt(u^2+v^2)
+  HProf GetProjResponse(int nbin=128, bool fgnorm1=true);
 
   double CumulResp(Four2DResponse& rd, double theta=0., double phi=0.);
@@ -143,5 +162,5 @@
   double AddToHisto(double kx0, double ky0, double x, double y, double w, bool fgfh);
 
-  double lambda_, dmax_;
+  double lambda_, dmax_, kmax_;
   vector<Dish> dishes_;
   bool fgnoauto_;
@@ -153,4 +172,5 @@
   //   Histo2D h2w_, h2cnt_;
   QHis2D h2w_;
+  bool fgcomputedone_;
   int mcnt_;
   int prtlev_,prtmodulo_;
Index: /trunk/Cosmo/RadioBeam/pknoise.cc
===================================================================
--- /trunk/Cosmo/RadioBeam/pknoise.cc	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/pknoise.cc	(revision 3947)
@@ -44,5 +44,6 @@
        << "    -renmax MaxValue (default : Do NOT renormalize 2D response value \n"   
        << "    -scut SCutValue (default=100.) \n"
-       << "    -ngen NGen (default=1) number of noise fourier amp generations \n"
+       << "    -ngen NGen (default=0) number of noise fourier amp generations \n"
+       << "       NGen=0 -> Call ComputeNoisePk(), else generate Fourier Amplitudes (random) \n"
        << "    -z redshift (default=0.7) \n"
        << "    -prt PrtLev,PrtModulo (default=0,10) " << endl;
@@ -74,5 +75,5 @@
   bool fgrenorm=false;
   double rmax=1.;
-  int NMAX = 1;
+  int NMAX = 0;
   double SCut=0.;
   double z_Redshift=0.7 ;  // 21 cm at z=0.7 -> 0.357 m  
@@ -141,7 +142,10 @@
       cout << "pknoise[1]: initializing Four2DRespTable from file" << resptblname << endl;
       resptbl.readFromPPF(resptblname);
+      cout << "pknoise[1.b] Four2DRespTable.LambdaRef=" << resptbl.getLambdaRef() << " setLambda(" 
+	   << lambda << ")" << endl;
+      resptbl.setLambda(lambda);
       arep_p=&resptbl;
       if (fgrenorm) {
-	cout << " pknoise[1.b] call to resptbl.renormalize(" << rmax << ")"; 
+	cout << "pknoise[1.c] call to resptbl.renormalize(" << rmax << ")"; 
 	double omax=resptbl.renormalize(rmax);
 	cout << " ... Old Max=" << omax << endl;
@@ -159,9 +163,25 @@
 
     cout << " pknoise[3]: Computing Noise P(k) using PkNoiseCalculator ..." << endl;
-    PkNoiseCalculator pkn(m3d, *(arep_p), SCut, NMAX, tits.c_str());
-    pkn.SetPrtLevel(prtlev,prtmod);
-    HProf hpn = pkn.Compute();
-    po << hpn;
- 
+    if (NMAX>0) {
+      PkNoiseCalculator pkn(m3d, *(arep_p), SCut, NMAX, tits.c_str());
+      pkn.SetPrtLevel(prtlev,prtmod);
+      HProf hpn = pkn.Compute();
+      po << hpn;
+    }
+    else {
+      Histo fracmodok;
+      DataTable dtnoise;
+      HProf hpn = m3d.ComputeNoisePk(*(arep_p),fracmodok,dtnoise,SCut);
+      po << hpn;
+      string outfile2 = "x"+outfile;
+      POutPersist po2(outfile2);
+      HProf h1dnoise=arep_p->GetProjNoiseLevel();
+      HProf h1drep=arep_p->GetProjResponse();
+      po2 << PPFNameTag("dtnoise") << dtnoise;
+      po2 << PPFNameTag("hpnoise") << hpn;
+      po2 << PPFNameTag("fracmodok") << fracmodok;
+      po2 << PPFNameTag("h1dnoise") << h1dnoise;
+      po2 << PPFNameTag("h1drep") << h1drep;
+    }
     rc = 0;
   }  // End of try bloc 
Index: /trunk/Cosmo/RadioBeam/plpkn.pic
===================================================================
--- /trunk/Cosmo/RadioBeam/plpkn.pic	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/plpkn.pic	(revision 3947)
@@ -17,4 +17,5 @@
 setup5 
 scalewz 0.7 2500 100 1
+# scalewz 0.25 989.13 80.26 1
 y1 = 500*$cct21
 y2 = 9e4*$cct21
Index: /trunk/Cosmo/RadioBeam/qhist.h
===================================================================
--- /trunk/Cosmo/RadioBeam/qhist.h	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/qhist.h	(revision 3947)
@@ -43,8 +43,20 @@
   double Add(r_8 x, r_8 y, r_8 w, bool fgfh);
   void SetZeroBin(r_8 x=0, r_8 y=0);
-  inline double WBinX() { return dxb; }
-  inline double WBinY() { return dyb; }
+  inline r_8 WBinX() { return dxb; }
+  inline r_8 WBinY() { return dyb; }
+  inline sa_size_t NBinX() { return nx; }
+  inline sa_size_t NBinY() { return ny; }
+  inline r_8 X(sa_size_t i) { return (dxb*((double)i+0.5)+xmin); }
+  inline r_8 Y(sa_size_t j) { return (dyb*((double)j+0.5)+ymin); }
+  inline r_8 operator()(sa_size_t i, sa_size_t j) { return aw(i,j); }
+
+  inline void GetMinMax(r_8& min, r_8& max) { return aw.MinMax(min,max); }
+
+  inline double SumWBinZero() { return sumw0; }
+  void Renormalize(r_8 mfac=1.) { aw*=mfac; }
+
   Histo2D Convert();
   
+protected:
   r_8 xmin,xmax,ymin,ymax;
   r_8 dxb,dyb;
Index: /trunk/Cosmo/RadioBeam/radutil.cc
===================================================================
--- /trunk/Cosmo/RadioBeam/radutil.cc	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/radutil.cc	(revision 3947)
@@ -45,5 +45,5 @@
   double zz = 1.+z_;
   double cc=zz*zz/sqrt(OmegaMatter_*zz*zz*zz+OmegaLambda_);
-  cc *= ((h100_/0.7)*(OmegaBaryons_/0.04)*(fracHI_/0.1));
-  return (cc*0.57);
+  cc *= ((h100_/0.7)*(OmegaBaryons_/0.044)*(fracHI_/0.01));
+  return (cc*0.054);
 }
Index: /trunk/Cosmo/RadioBeam/radutil.h
===================================================================
--- /trunk/Cosmo/RadioBeam/radutil.h	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/radutil.h	(revision 3947)
@@ -37,5 +37,5 @@
   // Definition des parametres cosmologiques utiles pour le calcul de la temperature d'emission a 21 cm
   // retourne la valeur de OmegaLambda (univers plat)
-  double setCosmoParam(double omegamatter=0.02581, double omegabaryon=0.0441, double h100=0.719, double fracHI=0.02);
+  double setCosmoParam(double omegamatter=0.2581, double omegabaryon=0.0441, double h100=0.719, double fracHI=0.02);
   inline void setFracHI(double fracHI=0.02) { fracHI_=fracHI; }
 
Index: /trunk/Cosmo/RadioBeam/repicon.cc
===================================================================
--- /trunk/Cosmo/RadioBeam/repicon.cc	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/repicon.cc	(revision 3947)
@@ -32,8 +32,8 @@
 
 // pour sauver la reponse mdresp et la config des dishes dans un fichier PPF 
-void SaveDTVecDishesH2Resp(string& outfile, vector<Dish>& vdishes, Four2DRespTable& mdresp);
+void SaveDTVecDishesH2Resp(POutPersist& po, vector<Dish>& vdishes, Four2DRespTable& mdresp);
 
 // ---------------------------------------------------------------------
-// main program for computing interferometer response un (u,v) plane
+// main program for computing interferometer response in (u,v) plane
 // R. Ansari  - Avril-Juin 2010 
 // ---------------------------------------------------------------------
@@ -42,18 +42,24 @@
 {
   cout << " Usage: repicon [-parname Value] configId OutPPFName \n" 
-       << " configIds: f4x4,f8x8,f20x20, confA,confB,confC,confD, hex12,cross11, f3cyl,f8cyl,f3cylp,f8cylp \n"
-       << "   f4x4 , f8x8 , f20x20 Filled array of nxn dishes  \n" 
+       << " configIds: f4x4,f8x8,f11x11,f20x20, confA,confB,confC,confD, hex12,cross11, \n" 
+       <<  "           f3cyl,f8cyl,f3cylp,f8cylp, nan12,nan24,nan36,nan40,nan128 \n"
+       << "   f4x4 , f8x8 , f11x11 , f20x20 Filled array of nxn dishes  \n" 
        << "   confA , confB, confC, confD : semi-filled array of dishes \n"
        << "   hex12,cross11 : ASKAP like double hexagonal (12xD=12m), cross config (11xD=12m) \n"
+       << "   nan12,nan24,nan36,nan40,nan128: 3,4,6,5,8 cylinder like configurations with 12,24,36,40,128 dishes\n"
        << "   f3cyl, f8cyl , f3cylp, f8cylp : filled array of non perfect/perfect of n cylinders \n"
        << "   f4cylw : filled array of 4 perfect of wide cylinders \n"
-       << "  [ -parname value] : -renmax -z -prt -D -lmax \n" 
+       << "   pit2cyl,pit2cylw : 2 cylinders 15mx7m, 20mx10 for z=0.3 at Pittsburg (perfect) \n"
+       << "   pit4cyl : 4 cylinders 32mx8m for z=0.3-0.7 at Pittsburg (perfect) \n"
+       << "  [ -parname value] : -renmax -z -prt -D -lmax -eta -autocor \n" 
        << "    -renmax MaxValue (default : Do NOT renormalize 2D response value \n"   
        << "    -z redshift (default=0.7) --> determines Lambda \n"
        << "    -D DishDiameter (default=5 m) \n"
+       << "    -eta fill_factor (default=0.90) \n"
        << "    -lmax array extension (default=100 m ) for response calculation kmax \n" 
        << "    -repnxy Nx,Ny (default=200,200) ResponseTable binning \n"
        << "    -beamnxy Nx,Ny (default=200,200) Beam sampling \n"
        << "    -rot ThetaMaxDeg,NTheta,PhiMaxDeg,NPhi (default NO-rotate/pointing -> 23,10,30,15 ) \n" 
+       << "    -autocor : keep antenna auto-correlation signals (default NO) \n"
        << "    -prt PrtLev,PrtModulo (default=0,10) \n" 
        << endl; 
@@ -82,4 +88,5 @@
   
   double Ddish=5.;
+  double Eta=0.9;
   bool fgDfixed=false;
   double LMAX = 100.;  // taille de la zone, pour calcul kmax 
@@ -94,4 +101,5 @@
   int NBX=200;
   int NBY=200;
+  bool fgnoauto=true;  // do not keep antenna autocorrelation signal 
   
   int ka=1;
@@ -106,4 +114,7 @@
       Ddish=atof(arg[ka+1]);  fgDfixed=true;  ka+=2;
     }
+    else if (strcmp(arg[ka],"-eta")==0) {
+      Eta=atof(arg[ka+1]);  ka+=2;
+    }
     else if (strcmp(arg[ka],"-lmax")==0) {
       LMAX=atof(arg[ka+1]);  fgLMAXfixed=true;  ka+=2;
@@ -117,4 +128,7 @@
     else if (strcmp(arg[ka],"-rot")==0) {
       sscanf(arg[ka+1],"%lg,%d,%lg,%d",&thetamxdeg,&nteta,&phimxdeg,&nphi);  fgpoint=true;  ka+=2;
+    }
+    else if (strcmp(arg[ka],"-autocor")==0) {
+      fgnoauto=false;   ka+=1;
     }
     else if (strcmp(arg[ka],"-prt")==0) {
@@ -143,5 +157,4 @@
     LAMBDA = conv.getLambda();
     
-    double Eta=0.95;
     double cylW=12.;   // Largeur des cylindres 
     double cylRL=2*LAMBDA;  // Longeur des elements de reception le long du cylindre
@@ -150,5 +163,6 @@
 
     vector<Dish> vdishes;
-    cout << " repicon[1] : creating dish vector/ MultiDish for configuration: " << config << endl;
+    cout << " repicon[1] : creating dish vector/ MultiDish for configuration: " << config 
+	 << " Ddish=" << Ddish << " m. Eta=" << Eta << endl;
 
     if (config=="f4x4") {
@@ -157,4 +171,7 @@
     else if (config=="f8x8") {
       vdishes=CreateFilledSqConfig(8,Ddish, Eta);
+    }
+    else if (config=="f11x11") {
+      vdishes=CreateFilledSqConfig(11,Ddish, Eta);
     }
     else if (config=="f20x20") {
@@ -172,4 +189,19 @@
       cylW=25.;   cylRL=3*LAMBDA; 
       vdishes=CreateFilledCylConfig(4, 100, cylW, cylRL, etaW, etaRL, false);
+    }
+    else if (config=="pit2cyl") {
+      etaW=0.9; etaRL=0.9;
+      cylW=7.;   cylRL=2*LAMBDA; 
+      vdishes=CreateFilledCylConfig(2, 32, cylW, cylRL, etaW, etaRL, false);
+    }
+    else if (config=="pit2cylw") {
+      etaW=0.9; etaRL=0.9;
+      cylW=9.;   cylRL=2*LAMBDA; 
+      vdishes=CreateFilledCylConfig(2, 40, cylW, cylRL, etaW, etaRL, false);
+    }
+    else if (config=="pit4cyl") {
+      etaW=0.9; etaRL=0.9;
+      cylW=8.;   cylRL=2*LAMBDA; 
+      vdishes=CreateFilledCylConfig(4, 64, cylW, cylRL, etaW, etaRL, false);
     }
     else if (config=="f8cyl") {
@@ -206,4 +238,19 @@
       vdishes=CreateDoubleHexagonConfig(Ddish);
     }
+    else if (config=="nan12") {
+      vdishes=CreateConfigNancay12(Ddish, Eta);
+    }
+    else if (config=="nan24") {
+      vdishes=CreateConfigNancay24(Ddish, Eta);
+    }
+    else if (config=="nan36") {
+      vdishes=CreateConfigNancay36(Ddish, Eta);
+    }
+    else if (config=="nan40") {
+      vdishes=CreateConfigNancay40(Ddish, Eta);
+    }
+    else if (config=="nan128") {
+      vdishes=CreateConfigNancay128(Ddish, Eta);
+    }
     else {
       cout << " NON valid configuration option -> exit" << endl;
@@ -212,8 +259,7 @@
         
     double Dol = LMAX/LAMBDA;
-    bool fgnoauto = true;
     
     cout << " repicon[1] : Lambda=" << LAMBDA << " LMAX= " << LMAX << " NDishes=" << vdishes.size() 
-	 << " D-Dish=" << Ddish << " m." << endl;
+	 << " D-Dish=" << Ddish << " m. " << ((fgnoauto)?" NO-":"With-") << "AutoCorrelation" << endl;
  
     MultiDish mdish(LAMBDA, LMAX, vdishes, fgnoauto);
@@ -230,4 +276,6 @@
     Histo2D hrep = mdish.GetResponse();
     PrtTim("Apres mdish.GetResponse()");
+    HProf h1dnoise = mdish.GetProjNoiseLevel();
+    HProf h1drep = mdish.GetProjResponse();
 
     Four2DRespTable mdresp(hrep, Dol, LAMBDA);
@@ -241,6 +289,9 @@
 
     string outfile2 = "hdt_"+outfile;
-    cout << " repicon[4] : saving H2D-response, multidish config to PPF file " << outfile2 << endl;
-    SaveDTVecDishesH2Resp(outfile2, vdishes, mdresp);
+    POutPersist po2(outfile2);
+    cout << " repicon[4] : saving H1D/H2D-response, multidish config to PPF file " << outfile2 << endl;
+    po2 << PPFNameTag("h1dnoise") << h1dnoise;
+    po2 << PPFNameTag("h1drep") << h1drep;
+    SaveDTVecDishesH2Resp(po2, vdishes, mdresp);
 
     rc = 0;
@@ -265,5 +316,5 @@
 
 /*-- Nouvelle-Fonction --*/
-void SaveDTVecDishesH2Resp(string& outfile, vector<Dish>& vdishes, Four2DRespTable& mdresp)
+void SaveDTVecDishesH2Resp(POutPersist& po, vector<Dish>& vdishes, Four2DRespTable& mdresp)
 {
   char* names[5]={"did","posx","posy","diam","diamy"};
@@ -285,5 +336,4 @@
   }
   Histo2D h2rep=mdresp.GetResponse();
-  POutPersist po(outfile);
   po << PPFNameTag("mdish") << ntvd;
   po << PPFNameTag("h2rep") << h2rep;
Index: /trunk/Cosmo/RadioBeam/specpk.cc
===================================================================
--- /trunk/Cosmo/RadioBeam/specpk.cc	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/specpk.cc	(revision 3947)
@@ -266,4 +266,78 @@
 }
 
+// Compute noise power spectrum as a function of wave number k 
+// No random generation, computes profile of noise sigma^2
+// cells with amp^2=re^2+im^2>s2cut are ignored
+// Output : noise power spectrum (profile histogram)
+HProf Four3DPk::ComputeNoisePk(Four2DResponse& resp, Histo& fracmodok, DataTable& dt, 
+			       double s2cut, int nbin, double kmin, double kmax)
+{
+  // The second half of the array along Y (matrix rows) contain
+  // negative frequencies
+  //  int nbh = sqrt(fourAmp.SizeX()*fourAmp.SizeX()+fourAmp.SizeY()*fourAmp.SizeY()/4.+fourAmp.SizeZ()*fourAmp.SizeY()/4.);
+  // The profile histogram will contain the mean value of noise sigma^2
+  // as a function of wave-number k = sqrt((double)(kx*kx+ky*ky))
+  //  if (nbin < 1) nbin = nbh/2;
+  if ((kmax<0.)||(kmax<kmin)) {
+    kmin=0.;
+    double maxx=fourAmp.SizeX()*dkx_;
+    double maxy=fourAmp.SizeY()*dky_/2;
+    double maxz=fourAmp.SizeZ()*dkz_/2;
+    kmax=sqrt(maxx*maxx+maxy*maxy+maxz*maxz);
+  }
+  if (nbin<2) nbin=128; 
+  HProf hp(kmin, kmax, nbin);
+  hp.SetErrOpt(false);
+  Histo hmcnt(kmin, kmax, nbin);
+  Histo hmcntok(kmin, kmax, nbin);
+
+  uint_8 nmodeok=0;
+  // fourAmp represent 3-D fourier transform of a real input array. 
+  // The second half of the array along Y and Z contain negative frequencies
+  double kxx, kyy, kzz;
+  // sa_size_t is large integer type  
+  // We ignore 0th term in all frequency directions ...
+  for(sa_size_t kz=1; kz<fourAmp.SizeZ(); kz++) {
+    kzz =  (kz > fourAmp.SizeZ()/2) ? (double)(fourAmp.SizeZ()-kz)*dkz_ : (double)kz*dkz_; 
+    for(sa_size_t ky=1; ky<fourAmp.SizeY(); ky++) {
+      kyy =  (ky > fourAmp.SizeY()/2) ? (double)(fourAmp.SizeY()-ky)*dky_ : (double)ky*dky_; 
+      for(sa_size_t kx=1; kx<fourAmp.SizeX(); kx++) {  // ignore the 0th coefficient (constant term)
+	double kxx=(double)kx*dkx_;
+	double wk = sqrt(kxx*kxx+kyy*kyy+kzz*kzz);
+	double amp2 = 1./resp(kxx, kyy); 
+	hmcnt.Add(wk);
+	if ((s2cut>1.e-9)&&(amp2>s2cut))  continue;
+	hmcntok.Add(wk);
+	hp.Add(wk, amp2);
+	nmodeok++;
+      }
+    }
+  }
+  if ((prtlev_>1)||((prtlev_>0)&&(s2cut>1.e-9))) {
+    cout << " Four3DPk::ComputeNoisePk/Info : NModeOK=" << nmodeok << " / NMode=" << fourAmp.Size() 
+	 << " -> " << 100.*(double)nmodeok/(double)fourAmp.Size() << "%" << endl;
+  }
+
+  fracmodok=hmcntok/hmcnt;
+
+  char* nomcol[5] = {"k","pnoise","nmode","nmodok","fracmodok"};
+  dt.Clear();
+  dt.AddDoubleColumn(nomcol[0]);    
+  dt.AddDoubleColumn(nomcol[1]);    
+  dt.AddIntegerColumn(nomcol[2]);    
+  dt.AddIntegerColumn(nomcol[3]);    
+  dt.AddFloatColumn(nomcol[4]);    
+  DataTableRow 	dtr = dt.EmptyRow();
+  for(int_4 ib=0; ib<hp.NBins(); ib++) {
+    dtr[0]=hp.BinCenter(ib);
+    dtr[1]=hp(ib);
+    dtr[2]=hmcnt(ib);
+    dtr[3]=hmcntok(ib);
+    dtr[4]=fracmodok(ib);
+    dt.AddRow(dtr);
+  }
+  return hp;
+}
+
 //-----------------------------------------------------
 // -- MassDist2D class :  2D mass distribution 
Index: /trunk/Cosmo/RadioBeam/specpk.h
===================================================================
--- /trunk/Cosmo/RadioBeam/specpk.h	(revision 3946)
+++ /trunk/Cosmo/RadioBeam/specpk.h	(revision 3947)
@@ -64,4 +64,7 @@
   HProf ComputePk(double s2cut=0., int nbin=256, double kmin=0., double kmax=-1.);
   void  ComputePkCumul(HProf& hp, double s2cut=0.);
+
+  HProf ComputeNoisePk(Four2DResponse& resp, Histo& fracmodok, DataTable& dt, 
+		       double s2cut=0., int nbin=256, double kmin=0., double kmax=-1.);
 
 protected:
