src/Pulses/PulGauss.cc
// This is a replica of the previous function that is used for debugging
// It contains additional print statements which indicate how the propagator
// is built step by step.  It isn't documented and may one day disappear.
This is a replica of the previous function that is used for debugging
It contains additional print statements which indicate how the propagator
is built step by step.  It isn't documented and may one day disappear.
This is a replica of the previous function that is used for debugging It contains additional print statements which indicate how the propagator is built step by step . It isn ' t documented and may one day disappear .
gen_op Gpulse_UX(gen_op& H0rot, gen_op& Fxy,
                          int N, double ang, double tp, double fact)
   {
   int debug=1;
   double dN = double(N);
   double tdiv = tp/dN;				
   double gamB1 = GgamB1(ang, tp, N, fact);	
   row_vector Gs = Gvect(gamB1, N, fact);	
   Fxy.Op_base(H0rot);				
   double thang, tang = 0;
   double gB1 = Gs.getRe(0);			
   gen_op H = H0rot;
   H += gB1*Fxy;				
   if(debug)
     {
     cout << "\n\t\tPulse Length = " << tp;
     cout << "\n\t\tStep Length = " << tdiv;
     cout << "\n\t\tSteps = " << N;
     cout << "\n\t\tStep 1 at " << gB1 << " Hz";
     cout << ", Angle " << gB1*tdiv*360. << " Degrees";
     tang += gB1*tdiv*360.;
     cout << ", Total Angle " << tang << " Degrees";
     }
   gen_op Utlow = prop(H, tdiv);		
   gen_op Uthigh = Utlow;			
   gen_op U, Ustep;				
   double tpx = tdiv;				
   int nstps = 1;				
   complex gb = gB1/dN;				
   int i;
   for(i=1; i<N-1-i; i++)			
     { 						
     gB1 = Gs.getRe(i);				
     H = H0rot;					
     H += gB1*Fxy;				
     Ustep = prop(H, tdiv);			
     Utlow &= Ustep;				
     Uthigh *= Ustep;				
     if(debug)
       {
       tpx += tdiv;				
       nstps++;					
       gb += Gs.get(i)/dN;			
       cout << "\n\t\tStep " << i+1 << " at " << gB1 << " Hz";
       cout << ", Angle " << gB1*tdiv*360. << " Degrees";
       tang += gB1*tdiv*360.;
       cout << ", Total Angle " << tang << " Degrees";
       }
     }
   thang = tang;				
   if(N == 1) U = Utlow;			
   else if(i == N-1-i)				
     {						
     gB1 = Gs.getRe(i);
     H = H0rot;
     H += gB1*Fxy;			
     U = prop(H,tdiv);
     U *= Utlow;
     U &= Uthigh;
     if(debug)
       {
       cout << "\n\t\tStep " << i+1 << " at " << gB1 << " Hz";
       cout << ", Angle " << gB1*tdiv*360. << " Degrees";
       tang += gB1*tdiv*360.;
       cout << ", Total Angle " << tang << " Degrees (MIDDLE)";
       cout << "\n\t\tEnd Steps Add an Additional Angle of " << thang;
       cout << " Degrees, Total Angle " << thang+tang << " Degrees";
       }
     }
   else
     {
     U = Uthigh*Utlow;			
     if(debug)
       {
       cout << "\n\t\tEnd Steps Add an Additional Angle of " << thang;
       cout << " Degrees, Total Angle " << thang+tang << " Degrees";
       }
     }
  if(debug)
    {
    if(N==1)
      {
      tpx = tdiv; 
      nstps = 1;
      gb =  Gs.get(0)/dN;
      }
    else if(i == N-1-i)
      {
      tpx = 2*tpx+tdiv; 
      nstps = 2*nstps+1;
      gb =  2*gb + Gs.get(i)/dN;
      }
    else
      {
      tpx *=2;
      nstps *= 2;
      gb *= 2.0;
      }
    cout << "\n\n\tSUPPOSED PULSE LENGTH = " << tp;
    cout << "\n\tSUMMED PULSE LENGTH = " << tpx;
    cout << "\n\tREQUIRED GAMMA B1 = " << gamB1;
    cout << "\n\tINTEGRATED GAMMA B1 = " << gb;
    cout << "\n\tSUPPOSED NUMBER OF STEPS = " << N;
    cout << "\n\tACTUAL NUMBER OF STEPS = " << nstps;
    cout << "\n\tSUPPOSED PULSE ANGLE = " << ang;
    cout << "\n\tACTUAL PULSE ANGLE = " << tpx*gb*360.0 << "\n";
    }
   return U;
   }
gen_op Gpulse_UX ( gen_op & H0rot , gen_op & Fxy , int N , double ang , double tp , double fact ) { int debug = 1 ; double dN = double ( N ) ; double tdiv = tp / dN ; double gamB1 = GgamB1 ( ang , tp , N , fact ) ; row_vector Gs = Gvect ( gamB1 , N , fact ) ; Fxy . Op_base ( H0rot ) ; double thang , tang = 0 ; double gB1 = Gs . getRe ( 0 ) ; gen_op H = H0rot ; H += gB1 * Fxy ; if ( debug ) { cout << " \n \t \t " << tp ; cout << " \n \t \t " << tdiv ; cout << " \n \t \t " << N ; cout << " \n \t \t " << gB1 << " " ; cout << " " << gB1 * tdiv * 360. << " " ; tang += gB1 * tdiv * 360. ; cout << " " << tang << " " ; } gen_op Utlow = prop ( H , tdiv ) ; gen_op Uthigh = Utlow ; gen_op U , Ustep ; double tpx = tdiv ; int nstps = 1 ; complex gb = gB1 / dN ; int i ; for ( i = 1 ; i < N - 1 - i ; i ++ ) { gB1 = Gs . getRe ( i ) ; H = H0rot ; H += gB1 * Fxy ; Ustep = prop ( H , tdiv ) ; Utlow &= Ustep ; Uthigh *= Ustep ; if ( debug ) { tpx += tdiv ; nstps ++ ; gb += Gs . get ( i ) / dN ; cout << " \n \t \t " << i + 1 << " " << gB1 << " " ; cout << " " << gB1 * tdiv * 360. << " " ; tang += gB1 * tdiv * 360. ; cout << " " << tang << " " ; } } thang = tang ; if ( N == 1 ) U = Utlow ; else if ( i == N - 1 - i ) { gB1 = Gs . getRe ( i ) ; H = H0rot ; H += gB1 * Fxy ; U = prop ( H , tdiv ) ; U *= Utlow ; U &= Uthigh ; if ( debug ) { cout << " \n \t \t " << i + 1 << " " << gB1 << " " ; cout << " " << gB1 * tdiv * 360. << " " ; tang += gB1 * tdiv * 360. ; cout << " " << tang << " " ; cout << " \n \t \t " << thang ; cout << " " << thang + tang << " " ; } } else { U = Uthigh * Utlow ; if ( debug ) { cout << " \n \t \t " << thang ; cout << " " << thang + tang << " " ; } } if ( debug ) { if ( N == 1 ) { tpx = tdiv ; nstps = 1 ; gb = Gs . get ( 0 ) / dN ; } else if ( i == N - 1 - i ) { tpx = 2 * tpx + tdiv ; nstps = 2 * nstps + 1 ; gb = 2 * gb + Gs . get ( i ) / dN ; } else { tpx *= 2 ; nstps *= 2 ; gb *= 2.0 ; } cout << " \n \n \t " << tp ; cout << " \n \t " << tpx ; cout << " \n \t " << gamB1 ; cout << " \n \t " << gb ; cout << " \n \t " << N ; cout << " \n \t " << nstps ; cout << " \n \t " << ang ; cout << " \n \t " << tpx * gb * 360.0 << " \n " ; } return U ; }
This is a replica of the previous function that is used for debugging
It contains additional print statements which indicate how the propagator
is built step by step.
This is a replica of the previous function that is used for debugging It contains additional print statements which indicate how the propagator is built step by step .
// Input 	H0rot   : Static Hamiltonian * //		Fxy	: RF-Field Hamiltonian * //              N       : Number of Gaussian steps //		ang     : Pulse rotation angle (degrees) //              tp      : Gaussian pulse length (sec) //		fact    : Cutoff factor [0.0,1.0] // Output       U       : Propagator for a Gaussian pulse // Note                 : HOrot is the static Hamiltonian //                        in the rotating frame of the pulse. // Note                 : Fxy implicitly contains rf-pulse phase //		U         = U       * U     * U //		 Gaussian    t>tp/2    tp/2    t<tp/2 // Inc. prop. time // Get field strength // RF intensities // Set Fxy basis for easy add // Field strength first step // Total Ham., 1st step // Propagator for 1st step // Propagator for Nth step // Working propagators // Used for debugging // Used for debugging // Used for debugging // Loop 1st half of Gaussian // steps, symmetry for 2nd half // RF-field strength this step // Static Hamiltonian // Total Hamiltonian, this step // Prop. for this pulse step // Utlow = Ustep*Utlow // Uthigh = Uthigh*Ustep // Used for debugging // Used for debugging // Used for debugging // Used for debugging // For one step, U is 1st step // For odd number steps (odd N) // need to still do middle step // Total Ham. , middle step // Even # steps (even N), this U