Index: trunk/Cosmo/SimLSS/cmvhshsns.cc
===================================================================
--- trunk/Cosmo/SimLSS/cmvhshsns.cc	(revision 3632)
+++ trunk/Cosmo/SimLSS/cmvhshsns.cc	(revision 3653)
@@ -1,4 +1,6 @@
 // lobes pour HSHS N-S
 // > cmvhshsns -n -d 0.105 -g 4,0.105,0. -i 2,0.41,0. -t 50,180,0 -f 1420 1410 1430
+// > cmvhshsns -n -d 0.105 -g 2,0.105,0. -i 2,0.41,0. -t 50,180,0 -f 1420 1410 1430
+// > cmvhshsns -n -D -1 -d 10 -g 1 -i 10,10,0. -t 50,180,0 -f 1420 1410 1430
 #include "sopnamsp.h"
 #include "machdefs.h"
@@ -20,8 +22,10 @@
 cout<<"cmvhshsns [...] val1 val2 ..."<<endl
   <<"  -f : val1... sont des frequences en MHz"<<endl
-  <<"       (defaut: longueurs d'onde en m"<<endl
+  <<"       (defaut: longueurs d'onde en m)"<<endl
   <<"  -n : hauteur des lobes normalises a 1"<<endl
   <<"  -d L : longueur totale dipole"<<endl
-  <<"         L==0. alors approximation du  dipole"<<endl
+  <<"         si \"-D -1\" largeur du cylindre (etude E-W)"<<endl
+  <<"  -D : +1 lobe du dipole approximation de Hertz"<<endl
+  <<"       -1 le lobe du dipole est remplace par le sinc() du cylindre (etude E-W)"<<endl
   <<"  -g N_g,D_g,Theta_g : regroupement de dipoles"<<endl
   <<"                       N_g<=1 pas de regroupement"<<endl
@@ -30,8 +34,8 @@
   <<"  -t Nang,Tmax,Tcent : nombre de pts entre 2 zeros consecutifs (def=25)"<<endl
   <<"                       angle maxi (deg, def=180), angle central (deg, def=0)"<<endl
-  <<"  -p : interprete Theta_{g,i} en picosecondes"<<endl
+  <<"  -p : interprete Theta_{g,i} en picosecondes de temps"<<endl
   <<"..distances L,D: >0 en m , <0 en unites de longeur d'onde"<<endl
   <<"..angles Theta: en deg"<<endl
-  <<"                en secondes si option \"-p\""<<endl
+  <<"                en picosecondes de temps si option \"-p\""<<endl
   <<endl;
 }
@@ -39,4 +43,5 @@
 double thetafromdt(double &theta,double dt,double lambda,double dconsec);
 double dtfromtheta(double theta,double lambda,double dconsec);
+double princpiclarge(double lambda,int N, double d,double sth=0.);
 
 //----------------------------------------------------------------
@@ -45,5 +50,5 @@
  const double torad = M_PI/180.;
 
- // --- longueur d'onde em m
+ // --- longueur d'onde en m
  vector<double> Lambda, Nu;
  bool argfreq = false;
@@ -53,4 +58,5 @@
  // --- dipole de longeur totale L (2 brins de L/2)
  double L_d = -0.5;   // >0 en m, <0 en unite de lambda
+ int dipolType = 0; // 1=dipole de Hertz, -1=sinc(), sinon dipole a 2 brins
 
  // --- groupes: regroupement des dipoles
@@ -79,5 +85,5 @@
  // Decodage des arguments
  char c;
- while((c = getopt(narg,arg,"hpnfd:g:i:t:")) != -1) {
+ while((c = getopt(narg,arg,"hpnfD:d:g:i:t:")) != -1) {
   switch (c) {
   case 'f' :
@@ -89,4 +95,7 @@
   case 'n' :
     normone = true;
+    break;
+  case 'D' :
+    dipolType = atoi(optarg);
     break;
   case 'd' :
@@ -127,5 +136,5 @@
  for(unsigned short i=0;i<Lambda.size();i++)printf(" %.3f m ,  %.3f MHz\n",Lambda[i],Nu[i]/1.e6);
 
- cout<<"Dipole : longueur totale L="<<L_d<<endl;
+ cout<<"Dipole : type="<<dipolType<<", longueur (ou largeur) totale L="<<L_d<<endl;
  if(L_d==0.) return -4;
 
@@ -158,5 +167,5 @@
    cout<<"\n\n>>> Lambda = "<<lambda<<" m ,  nu = "<<nu/1.e6<<" MHz"<<endl;
    double ld = (L_d<0.) ? -L_d*lambda : L_d;
-   cout<<"dipole: ld="<<ld<<" m"<<endl;
+   cout<<"dipole: ld="<<ld<<" m  (type="<<dipolType<<")"<<endl;
    double dg = (D_g<0.) ? -D_g*lambda : D_g;
    cout<<"groupe: ("<<N_g<<"), dg="<<dg<<" m -> "<<dg*N_g<<" m"<<endl;
@@ -181,5 +190,15 @@
    double sthi = sin(thi);
 
-   //... distance approx entre 2 zeros
+   //... largeur du pic principale (distance entre les 2 zeros de part et d'autre)
+   if(N_g>1) {
+     double t = princpiclarge(lambda,N_g,dg,sthg);
+     if(t>0.) cout<<"groupe: largeur du lobe principal: "<<t<<" rad   "<<t/torad<<" deg"<<endl;
+   }
+   if(N_i>1) {
+     double t = princpiclarge(lambda,N_i,di,sthi);
+     if(t>0.) cout<<"interf: largeur du lobe principal: "<<t<<" rad   "<<t/torad<<" deg"<<endl;
+   }
+
+   //... distance approx entre 2 zeros (hors pic principale)
    double dzero = M_PI/2.;
    if(N_g>1) {
@@ -204,4 +223,5 @@
    double dt = Tmax*torad/npt;
    cout<<"nombre de points dans la boucle "<<2*npt+1<<" , dt="<<dt/torad<<endl;
+   double vhmax = -1.e20, thmax=0.;
    for(int i=-npt;i<=npt;i++) {
      double t = Tcent*torad + i*dt;
@@ -210,8 +230,12 @@
      double deltag = M_PI*dg/lambda*(st-sthg);
      double deltai = M_PI*di/lambda*(st-sthi);
-     double ant = (ld==0.) ? AntDipole(ld/lambda,ta): AntCentFed(ld/lambda,ta);
+     double ant;
+     if(dipolType<0)        ant = LobeSinc(ld/lambda,t);
+       else if(dipolType>0) ant = AntDipole(ld/lambda,ta);
+         else               ant = AntCentFed(ld/lambda,ta);
      double intfg = (N_g==1) ? 1.: SinNXsX_Sqr(deltag,N_g)/norme_g;
      double intfi = (N_i==1) ? 1.: SinNXsX_Sqr(deltai,N_i)/norme_i;
      double intf = ant*intfg*intfi;
+     if(fabs(t)<M_PI/2. && intf>vhmax) {vhmax=intf; thmax=t;}
 
      xnt[0] = t/torad;
@@ -222,4 +246,5 @@
      nt.Fill(xnt);
    }
+   cout<<"Found maxi ("<<vhmax<<") for t = "<<thmax<<" rad = "<<thmax/torad<<" deg"<<endl;
 
    //...ecriture ppf
@@ -229,8 +254,9 @@
    DVList dvl;
    dvl("Lambda") = lambda; dvl("Nu") = nu;
-   dvl("Ld") = ld;
-   dvl("Ng") = N_g; dvl("Dg") = dg; dvl("Thg") = thg; dvl("Tg") = tg;
-   dvl("Ni") = N_i; dvl("Di") = di; dvl("Thi") = thi; dvl("Ti") = ti;
+   dvl("Ld") = ld; dvl("dipolType") = dipolType;
+   dvl("Ng") = N_g; dvl("Dg") = dg; dvl("Thg") = thg/torad; dvl("Tg") = tg;
+   dvl("Ni") = N_i; dvl("Di") = di; dvl("Thi") = thi/torad; dvl("Ti") = ti;
    dvl("Tmax") = Tmax; dvl("Tcent") = Tcent;
+   dvl("vhmax") = vhmax; dvl("thmax") = thmax/torad;
    sprintf(str,"dvl_%d",il);
    pos << PPFNameTag(str) << dvl;
@@ -238,37 +264,28 @@
 
    //... remplissage des zeros et des maximas principaux
-   {
-   const int nnt=2; double xnt[nnt];
-   const char *namev[nnt] = {"t","ztyp"};
-   if(N_g>1) {
+   for(int i=0;i<2;i++) {  // i==0 groupes , i==1 interferences
+     int N = N_g;
+     double d = dg, sth = sthg, norme = norme_g;
+     sprintf(str,"ntzg_%d",il);
+     if(i==1) {
+       N = N_i;
+       d = di; sth = sthi; norme = norme_i;
+       sprintf(str,"ntzi_%d",il);
+     }
+     if(N<=1) continue;
+     const int nnt=2; double xnt[nnt];
+     const char *namev[nnt] = {"t","ztyp"};
      NTuple nt(nnt,namev);
-     double no = (normone) ? 1.: norme_g;
-     for(int is=-1;is<=1;is+=2) {
+     double no = (normone) ? 1.: norme;
+     for(int is=-1;is<=1;is+=2) {  // ki>=0 ou k<0
        int k0 = (is==1) ? 0: -1;
        for(int k=k0;;k+=is) {
-         xnt[0] =  k*lambda/(N_g*dg)+sthg;
-         xnt[1] = (k%N_g==0) ? no: 0.;
+         xnt[0] =  k*lambda/(N*d)+sth;
+         xnt[1] = (k%N==0) ? no: 0.;  // zero ou maxi?
          if(fabs(xnt[0])>1.) break;
          xnt[0]=asin(xnt[0])/torad; nt.Fill(xnt);
        }
      }
-     sprintf(str,"ntzg_%d",il);
      pos << PPFNameTag(str) << nt;
-   }
-   if(N_i>1) {
-     NTuple nt(nnt,namev);
-     double no = (normone) ? 1.: norme_i;
-     for(int is=-1;is<=1;is+=2) {
-       int k0 = (is==1) ? 0: -1;
-       for(int k=k0;;k+=is) {
-         xnt[0] =  k*lambda/(N_i*di)+sthi;
-         xnt[1] = (k%N_i==0) ? no: 0.;
-         if(fabs(xnt[0])>1.) break;
-         xnt[0]=asin(xnt[0])/torad; nt.Fill(xnt);
-       }
-     }
-     sprintf(str,"ntzi_%d",il);
-     pos << PPFNameTag(str) << nt;
-   }
    }
 
@@ -331,5 +348,23 @@
 }
 
+double princpiclarge(double lambda,int N, double d,double sth)
+// Input:
+//   lambda : longueur d'onde
+//   N : nombre de dipole dans le regroupement
+//   d : distance entre deux dipoles
+//   sth : sin(theta) equivalent pour dephasage "electronique"
+// Return: si >0 distance angulaire en radians
+//         si -1 zero haut (droite) non-existant
+//         si -2 zero bas (gauche) non-existant
+{
+ double st1 = lambda/(N*d)+sth;
+ if(fabs(st1)>1.) return -1.;
+ double st2 = -lambda/(N*d)+sth;
+ if(fabs(st2)>1.) return -2.;
+ return fabs(asin(st1)-asin(st2));
+}
+
 /******************************
+delobjs *
 openppf cmvhshsns.ppf
 
@@ -342,5 +377,7 @@
 set cut -90<t&&t<90
 
-n/plot $nt.$t%_nl
+c++exec cout<<dvl_${l}<<endl;
+
+n/plot nt_$l.$t%_nl
 
 n/plot nt_$l.ant%$t $cut ! "nsta cpts green"
