Changeset 3653 in Sophya for trunk/Cosmo

Timestamp:

Jul 16, 2009, 12:33:57 PM (16 years ago)

Author:

cmv

Message:

add hshs EW, cmv 16/07/2009

Location:

trunk/Cosmo/SimLSS

Files:

• Makefile (modified) (3 diffs)
• cmvhshsew.cc (added)
• cmvhshsns.cc (modified) (18 diffs)
• geneutils.cc (modified) (2 diffs)
• geneutils.h (modified) (1 diff)

trunk/Cosmo/SimLSS/Makefile

@@ -21 +21 @@
        $(EXE)cmvtstsch $(EXE)cmvtstblack $(EXE)cmvtvarspec $(EXE)cmvdefsurv \
        $(EXE)cmvtintfun $(EXE)cmvconcherr $(EXE)cmvtinterp $(EXE)cmvtstagn \
-       $(EXE)cmvschdist $(EXE)cmvhshsns
+       $(EXE)cmvschdist $(EXE)cmvhshsns $(EXE)cmvhshsew
  
 PROGSOBJ = \
@@ -29 +29 @@
           $(OBJ)cmvtintfun.o $(OBJ)cmvtinterp.o \
           $(OBJ)cmvconcherr.o $(OBJ)cmvtluc.o $(OBJ)cmvtstagn.o $(OBJ)cmvschdist.o \
-          $(OBJ)cmvhshsns.o
+          $(OBJ)cmvhshsns.o $(OBJ)cmvhshsew.o
  
 LIBROBJ = \
@@ -280 +280 @@
         $(CXXCOMPILE) $(CXXREP) -I$(MYEXTINC) -o $@ cmvhshsns.cc
 
-##############################################################################
+cmvhshsew: $(EXE)cmvhshsew
+        echo $@ " done"
+$(EXE)cmvhshsew: $(OBJ)cmvhshsew.o $(LIBR)
+        $(CXXLINK) $(CXXREP) -o $@ $(OBJ)cmvhshsew.o $(MYLIB) 
+$(OBJ)cmvhshsew.o: cmvhshsew.cc
+        $(CXXCOMPILE) $(CXXREP) -I$(MYEXTINC) -o $@ cmvhshsew.cc
+
+##############################################################################

trunk/Cosmo/SimLSS/cmvhshsns.cc

@@ -1 +1 @@
 // 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 +22 @@
 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 +34 @@
   <<"  -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 +43 @@
 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 +50 @@
  const double torad = M_PI/180.;
 
- // --- longueur d'onde em m
+ // --- longueur d'onde en m
  vector<double> Lambda, Nu;
  bool argfreq = false;
@@ -53 +58 @@
  // --- 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 +85 @@
  // 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 +95 @@
   case 'n' :
     normone = true;
+    break;
+  case 'D' :
+    dipolType = atoi(optarg);
     break;
   case 'd' :
@@ -127 +136 @@
  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 +167 @@
    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 +190 @@
    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 +223 @@
    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 +230 @@
      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 +246 @@
      nt.Fill(xnt);
    }
+   cout<<"Found maxi ("<<vhmax<<") for t = "<<thmax<<" rad = "<<thmax/torad<<" deg"<<endl;
 
    //...ecriture ppf
@@ -229 +254 @@
    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 +264 @@
 
    //... 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 +348 @@
 }
 
+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 +377 @@
 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"

trunk/Cosmo/SimLSS/geneutils.cc

@@ -398 +398 @@
          d=8cm, distance entre 2 centres consecutifs 10cm
 #
-# ---- AntCentFed et AntDipole
+# ---- AntCentFed , AntDipole et LobeSinc
 # Input:
 #   t = angle entre la direction de la radiation et le fil de l'antenne (radian)
@@ -432 +432 @@
  double vd = pil*pil*st/2.;  vd *= vd;
  return vd;
+}
+
+double LobeSinc(double L, double trad)
+{
+ double x = M_PI*L*sin(trad);
+ return SinXsX_Sqr(x);
 }
 

trunk/Cosmo/SimLSS/geneutils.h

@@ -96 +96 @@
 double AntCentFed(double L, double trad);
 double AntDipole(double L, double trad);
+double LobeSinc(double L, double trad);
 
 //----------------------------------------------------