/* fmcpmghsqc.c

    This pulse sequence will allow one to perform the following experiment:

    2D CPMG-HSQC  Heteronuclear Single Quantum Coherence with enhanced sensitivity
                  for exchange-broadend signals.
                  Sequence employs PFGSE and dephasing of water.
 
    This sequence uses the standard three channel configuration
         1)  1H             - carrier (tof) @ 7.9 ppm [centre of amides]
         2) 13C             - carrier (dof) @ 117 ppm [centre of CO and Ca]
         3) 15N             - carrier (dof2)@ 119 ppm [centre of amide 15N]
 
    Set dm = 'nnnn', dmm = 'cccc'
    Set dm2 = 'nnny', dmm2 = 'cccp' [15N decoupling during acquisition]
                      dseq2 = 'waltz16'
 
    Must set phase = 1,2 for States-TPPI acquisition in t1 [15N] with
    gradient selection of 15N magnetization. [The fids must be manipulated
    (add/subtract) with 'grad_sort_nd' program (or equivalent) prior to
    regular processing.]
 
    Flags
        fsat            'y' for presaturation of H2O
        fscuba          'y' to apply scuba pulse after presaturation of H2O
        mess_flg        'y' for Messerlie type purging pulse
        f1180           'y' for 180 deg linear phase correction in F1
                            otherwise 0 deg linear phase correction
        c180_flg        'y' for 15N and 13C labelled samples
 
        Standard Settings
        fsat='n',fscuba='n',mess_flg='n',f1180='n',c180_flg='n'
          
        This CPMG version written by RM, Dec 29, 97
        (based on hsqc_gd_sl_seduce_500.c L. E. Kay July 7 1992)

        Adaptation by FM, Apr 7,1999
        1) altered gradients for coherence selection to be same
           as fmcpmghsqc_ref (=hsqc_gd_sl_seduce_500 of April 1999)
        2) water is dephased at the start of the sequence by an
           off-resonance pulse

        This version ONLY splits encoding gradient during t1 for 15N
        coherence selection.

        REF: F.A.A. Mulder et al.
             J. Biomol. NMR (1996) ....

        REF: L. E. Kay, P. Keifer, and T. Saarinen
             J. Am. Chem. Soc. 114, 10663-10665 (1992).

        REF: Zhang et. al.  
             J. Biomol. NMR 4, 845 (1994) 

*/

#include <standard.h>

static double d2_init = 0.0;

static int phi1[1] = {0},
           phi2[1] = {1},
           phi3[4] = {0,1,2,3}, 
           phi4[1] = {0},
           phi5[1] = {0},
           phi7[1] = {0},
           rec[2] = {0,2};

cpmg()
{
  double pwx2,jxh,tauxh,tau16;
  
  pwx2=getval("pwx2");
  jxh = getval("jxh");
  tauxh = 1.0/(4.0*jxh);
  tau16 = tauxh/16.0 - pwx2;

  initval(0.0,v1);
  initval(5.0,v3);
  initval(4.0,v5);

  loop(v5,v6);
   mod4(v1,v2);
   mod4(v3,v4);
   delay(tau16);
   sim3pulse(2*pw,0.0,2*pwx2,v2,zero,v2,0.0,0.0);
   delay(tau16);
   delay(tau16);
   sim3pulse(2*pw,0.0,2*pwx2,v4,zero,v4,0.0,0.0);
   delay(tau16);
   delay(tau16);
   sim3pulse(2*pw,0.0,2*pwx2,v2,zero,v2,0.0,0.0);
   delay(tau16);
   delay(tau16);
   sim3pulse(2*pw,0.0,2*pwx2,v4,zero,v4,0.0,0.0);
   delay(tau16);
   incr(v1);
   decr(v3);
  endloop(v6);
}

pulsesequence()
{
/* DECLARE VARIABLES */

 char        fscuba[MAXSTR],f1180[MAXSTR],fsat[MAXSTR],mess_flg[MAXSTR],c180_flg[MAXSTR],shp_sl[MAXSTR];

 int         phase,icosel,t1_counter;

 double      hscuba,                /* length of 1/2 scuba delay */
             tauxh,                 /* 1 / 4J(XH)                */
             tau16,               
             pwx2,                  /* PW90 for N-nuc            */
             pwx1,                  /* PW90 for C-nuc            */
             tsatpwr,               /* low power level for presat*/
             tpwrmess,              /* power level for Messerlie purge */
             dly_pg,                /* Messerlie purge delay         */
             dhpwr2,                /* power level for N hard pulses */
             dhpwr,                 /* power level for C hard pulses */
             jxh,                   /* coupling for XH           */
             tau1,                  /* t1/2 */
             sw1,
             pw_sl,                 /* pw for selective pulse at twprsl  */
             tpwrsl,

             gzlvl0,               /* level of grad. */
             gt0,                  /* grad time */
             gzlvl1,               /* level of grad. */
             gt1,                  /* grad time */
             gzlvl2,
             gt2,
             gzlvl3,
             gt3,
             gzlvl4,
             gt4,
             gzlvl5,
             gt5,
             BigT,                  /* constant time for N evolution */
             BigT1;


/* LOAD VARIABLES */

  jxh = getval("jxh");
  pwx2 = getval("pwx2");
  pwx1 = getval("pwx1");
  tsatpwr = getval("tsatpwr");
  tpwrmess = getval("tpwrmess");
  dly_pg = getval("dly_pg");
  dhpwr2 = getval("dhpwr2"); 
  dhpwr = getval("dhpwr"); 
  hscuba = getval("hscuba");
  phase = (int) (getval("phase") + 0.5);
  sw1 = getval("sw1");
  BigT = getval("BigT");
  BigT1 = getval("BigT1");
  pw_sl = getval("pw_sl");
  tpwrsl = getval("tpwrsl");

  gt0 = getval("gt0");
  gt1 = getval("gt1");
  gt2 = getval("gt2");
  gt3 = getval("gt3");
  gt4 = getval("gt4");
  gt5 = getval("gt5");
  gzlvl0 = getval("gzlvl0");
  gzlvl1 = getval("gzlvl1");
  gzlvl2 = getval("gzlvl2");
  gzlvl3 = getval("gzlvl3");
  gzlvl4 = getval("gzlvl4");
  gzlvl5 = getval("gzlvl5");

  getstr("fscuba",fscuba); 
  getstr("fsat",fsat); 
  getstr("f1180",f1180); 
  getstr("mess_flg",mess_flg);
  getstr("c180_flg",c180_flg);
  getstr("shp_sl",shp_sl);

/* check validity of parameter range */

    if((dm[A] == 'y' || dm[B] == 'y' || dm[C] == 'y' || dm[D] == 'y'))
        {
        printf("incorrect Dec1 decoupler flags!  ");
        abort(1);
    } 

    if((dm2[A] == 'y' || dm2[B] == 'y' || dm2[C] == 'y' ))
        {
        printf("incorrect Dec2 decoupler flags!  ");
        abort(1);
    } 

    if( tsatpwr > 8 )
    {
        printf("tsatpwr too large !!!  ");
        abort(1);
    }

    if ( tpwrmess > 57 )
    {
       printf("tpwrmess too high !!!");
       abort(1);
    }

    if (dly_pg > 0.010 )
    {
       printf("dly_pg too long !! ");
       abort(1);
    }

    if( dpwr > 50 )
    {
        printf("don't fry the probe, dpwr too large!  ");
        abort(1);
    }

    if( dpwr2 > 50 )
    {
        printf("don't fry the probe, dpwr2 too large!  ");
        abort(1);
    }

    if(gt0 > 15.0e-3 || gt1 > 15.0e-3 || gt2 > 15.0e-3 || gt3 > 15.0e-3 || 
       gt4 > 15.0e-3 || gt5 > 15.0e-3)
    {
        printf("gti must be less than 15 ms \n");
        abort(1);
    }

    if(mess_flg[A] == 'y') {
        printf("mess_flg should REALLY be set to n for optimum performance\n");
    }


/* LOAD VARIABLES */

  settable(t1, 1, phi1);
  settable(t2, 1, phi2);
  settable(t3, 4, phi3);
  settable(t4, 1, phi4);
  settable(t5, 1, phi5);
  settable(t7, 1, phi7);
  settable(t6, 2, rec);

/* INITIALIZE VARIABLES */

  tauxh = ((jxh != 0.0) ? 1/(4*(jxh)) : 2.25e-3);
  tau16 = tauxh/16.0;

/* Phase incrementation for hypercomplex data */

   if ( phase == 1 )     /* Hypercomplex in t1 */
   {
        tsadd(t5, 2, 4); 
        icosel = 1; /*  and reverse the sign of the gradient  */
   } 
   else icosel = -1;


/* calculate modification to phases based on current t1 values
   to achieve States-TPPI acquisition */
 
 
   if(ix == 1)
      d2_init = d2;
      t1_counter = (int) ( (d2-d2_init)*sw1 + 0.5);

      if(t1_counter %2) {
        tsadd(t7,2,4);
        tsadd(t6,2,4);
      }


/* set up so that get (-90,180) phase corrects in F1 if f1180 flag is y */

   tau1 = d2;
   if(f1180[A] == 'y')  tau1 += ( 1.0/(2.0*sw1) );
   tau1 = tau1/2.0;
   

/* BEGIN ACTUAL PULSE SEQUENCE */

status(A);
   rlpower(tsatpwr,TODEV);            /* Set power for presaturation  */
   rlpower(dhpwr,DODEV);              /* Set decoupler1 power to dhpwr */
   rlpower(dhpwr2,DO2DEV);            /* Set decoupler2 power to dhpwr2 */

/* Presaturation Period */

status(B);

/* option for Messerlie purge */

  if(mess_flg[A] == 'y') {
    rlpower(tpwrmess,TODEV);
    rgpulse(dly_pg,zero,20.0e-6,20.0e-6);
    rgpulse(dly_pg/1.62,one,20.0e-6,20.0e-6);
    rlpower(tsatpwr,TODEV);
 }

 if(fsat[0] == 'y')
{
  txphase(zero);
  rgpulse(d1,zero,20.0e-6,20.0e-6);
  rlpower(tpwr,TODEV);                /* Set power for hard pulses  */

    if (fscuba[0] == 'y')            /* Scuba pulse sequence */
    {  
      hsdelay(hscuba);

      rgpulse(pw,zero,1.0e-6,0.0);      /* 90x180y90x */
      rgpulse(2*pw,one,1.0e-6,0.0);
      rgpulse(pw,zero,1.0e-6,0.0);
 
      txphase(zero);
      delay(hscuba);        
    }
}

else  {
 rlpower(tpwr,TODEV);                /* Set power for hard pulses  */
 delay(d1);
}

status(C);

  rcvroff();
  delay(20.0e-6);

  /* eliminate all magnetization originating on 15N */

  dec2rgpulse(pwx2,zero,0.0,0.0);
  delay(2.0e-6);
  rgradient('z',gzlvl0);
  delay(gt0);
  rgradient('z',0.0);
  delay(150.0e-6);

  /* shaped pulse applied off-resonance to dephase water */
  rlpower(tpwrsl,TODEV);
  shaped_pulse(shp_sl,pw_sl,zero,2.0e-6,0.0);
  delay(2.0e-6);
  rlpower(tpwr,TODEV);
  /* end of shaped pulse */
  delay(2.0e-6);
  rgradient('z',gzlvl0);
  delay(gt0);
  rgradient('z',0.0);
  delay(150.0e-6);

  rgpulse(pw,zero,1.0e-6,0.0);    

  txphase(zero);
  dec2phase(zero);

  cpmg();

  rgpulse(pw,t2,0.0,0.0);

  delay(2.0e-6);
  rgradient('z',gzlvl2);
  delay(gt2);
  rgradient('z',0.0);
  delay(150.0e-6);

  dec2rgpulse(pwx2,t7,0.0,0.0);

  txphase(t4); dec2phase(t3);

  delay(tau1);                  /* delay=t2/2      */ 

  if(c180_flg[A] == 'y')
     simpulse(2*pw,2*pwx1,t4,zero,0.0,0.0);
  else
     rgpulse(2*pw,t4,0.0,0.0);

  txphase(zero);

  delay(tau1);

  delay(8.0e-6);
  rgradient('z',-1.0*gzlvl3);
  delay(gt3);
  rgradient('z',0.0);
  delay(2.0e-6);

  delay(BigT - gt3 - 2.0*GRADIENT_DELAY - 10.0e-6);

  dec2rgpulse(2*pwx2,t3,0.0,0.0);

  dec2phase(t5);

  delay(2.0e-6);
  rgradient('z',gzlvl3);
  delay(gt3);
  rgradient('z',0.0);
  delay(2.0e-6);
 
  if( (c180_flg[A] == 'y') && (pwx1 > pw) )
     delay(BigT + 2*pwx1 - gt3 - 2*GRADIENT_DELAY - 4.0e-6);  
  else 
     delay(BigT + 2*pw - gt3 - 2*GRADIENT_DELAY - 4.0e-6);  

  sim3pulse(pw,0.0e-6,pwx2,zero,zero,t5,0.0,0.0);   /* X read pulse */

  txphase(zero);
  dec2phase(zero);

  cpmg();

  sim3pulse(pw,0.0e-6,pwx2,one,zero,one,0.0,0.0);

  dec2phase(zero); txphase(zero);

  cpmg();

  rgpulse(pw,zero,0.0,0.0);

  delay(BigT1);

  rgpulse(2*pw,zero,0.0,0.0);

  delay(2.0e-6);
  rgradient('z',icosel*gzlvl4);
  delay(gt4);
  rgradient('z',0.0);
  delay(2.0e-6);

  rlpower(dpwr2,DO2DEV); /* low power for decoupling on dec channel 2 */ 
  rlpower(dpwr,DODEV); /* turn power down on dec channel */

  delay(BigT1 - gt4 - 4.0e-6 - 2.0*GRADIENT_DELAY - 2.0*POWER_DELAY);

/* acquire data */

status(D);
     setreceiver(t6);
}