/* hmqcj_lek_500_v2.c

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

       2D hmqcj for the measurement of NH-HA coupling constants.
                 F1        15N
                 F2(acq)   NH
      
       Uses channels 1 and 3 of the standard 3 channel configuration
       1) 1H        - carrier (tof) @ 4.7ppm [H2O]
       2)
       3) 15N       - carrier (dof2) @ 119ppm [ centre of amide N]

       Set dm = 'nnnn'  dmm = 'cccc'
           dm2 = 'nnny' dmm2 = 'cccp' [15N decoupling during acquisition]
                        dseq2 ='waltz16'

       Must set phase=1,2 for States-TPPI acquisition in t1 [15N]

       Flags
            fsat      'y' for presaturation of H2O
            fscuba    'y' to apply scuba pulse after presaturation of H2O
            mess_flg  'y' for 1H (H2O) purging before the relaxation delay
            f1180     'y' for 180deg linear phase correction in F1
            shp_flg   'y' for selective excitation and dephasing H2O mag.

       standard settings:
               fsat='n', fscuba='n', mes_flg='n', f1180='y', shp_flg='y'
   
       Written by L. E. Kay  05-12-92
       Modified by L. E. Kay 07-16-93 to encorporate gradients and boban

       Modified by L. E. Kay Sept 12/93 to move the rcvroff
       Modified by L.E.K on Feb, 22/94 to put in water selective pulses
       Modified by L.E.K on July 29 to eliminate the need for presat

 REF: L. E. Kay and A. Bax     J. Magn. Reson. 86, 110-126 (1990).

BMRB Pulse Sequence Accession Number: 55
*/

#include <standard.h>
#define PI 3.1416

static double d2_init = 0.0;

static int phi1[1] = {0},
           phi2[2] = {0,2},
           phi3[8] = {0,0,1,1,2,2,3,3}, 
           phi4[16] = {0,0,0,0,0,0,0,0,2,2,2,2,2,2,2,2},
           rec[16] = {0,2,2,0,0,2,2,0,
                      2,0,0,2,2,0,0,2};


pulsesequence()
{
/* DECLARE VARIABLES */

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


 int	     phase,
             t1_counter;

 double      hscuba,                /* length of 1/2 scuba delay */
             tauxh,                 /* 1 / 4J(XH)                */
             pwx2,                  /* PW90 for X-nuc on channel 3  */
             tsatpwr,               /* low power level for presat */
             dhpwr2,                /* power level for X hard pulses */
             jxh,                   /* coupling for XH           */
	     tau1,	      	    /* t1/2 */
	     sw1,
             tpwrmess,    /* power level for Messerlie purge */
             dly_pg1,     /* duration of first part of purge */
             tpwrsl,
             pw_sl,
             d1_wait,
 
             gzlvl1,
             gt1,
             gzlvl2,
             gt2;


/* LOAD VARIABLES */

  jxh = getval("jxh");
  pwx2 = getval("pwx2");
  tsatpwr = getval("tsatpwr");
  dhpwr2 = getval("dhpwr2"); 
  hscuba = getval("hscuba");
  phase = (int) (getval("phase") + 0.5);
  sw1 = getval("sw1");
  tpwrmess = getval("tpwrmess");
  dly_pg1 = getval("dly_pg1");
  d1_wait = getval("d1_wait");
  tpwrsl = getval("tpwrsl");
  pw_sl = getval("pw_sl");

  gzlvl1 = getval("gzlvl1");
  gt1 = getval("gt1");
  gzlvl2 = getval("gzlvl2");
  gt2 = getval("gt2");
 

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

/* check validity of parameter range */

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

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

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

    if( dpwr > 10 )
    {
	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( dhpwr2 > 62 )
    {
	printf("don't fry the probe, dhpwr2 too large!  ");
	abort(1);
    }

    if(tpwrmess > 55)
    {
        printf("dont fry the probe, tpwrmess too large < 55");
        abort(1);
    }

    if(dly_pg1 > 0.010)
    {
        printf("dont fry the probe, dly_pg1 is too long < 10 ms");
        abort(1);
    }

    if(gt1 > 15.0e-3 || gt2 > 15.0e-3) 
    {
        printf("gti is too long\n");
        abort(1);
    }

    if(tpwrsl > 25) 
    {
        printf("tpwrsl is too high\n");
        abort(1);
    }
   
/* LOAD VARIABLES */

  settable(t1, 1, phi1);
  settable(t2, 2, phi2);
  settable(t3, 8, phi3);
  settable(t4, 16, phi4);
  settable(t5, 16, rec);

/* INITIALIZE VARIABLES */

  tauxh = ((jxh != 0.0) ? 1/(2*(jxh)) : 3.57e-3);

/* Phase incrementation for hypercomplex data */

   if ( phase == 2 )     /* Hypercomplex in t1 */
   {
        tsadd(t2, 1, 4);
   } 

/* calculate modifications 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(t2,2,4);
       tsadd(t5,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) - 2.0*pw - 4.0/PI*pwx2 );
   else tau1 = tau1 - 2.0*pw - 4.0/PI*pwx2;
   if(tau1 < 0.2e-6) tau1 = 0.4e-6;
   tau1 = tau1/2.0;
   

/* BEGIN ACTUAL PULSE SEQUENCE */

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

/* Presaturation Period */

status(B);

   if (mess_flg[A] == 'y') {

      rlpower(tpwrmess,TODEV);
      rgpulse(dly_pg1,zero,2.0e-6,2.0e-6);
      rgpulse(dly_pg1/1.62,one,2.0e-6,2.0e-6);

      rlpower(tsatpwr,TODEV);
   }

 if(fsat[0] == 'y')
{
  txphase(zero);
  
  if(d1_wait == 0.0)
    rgpulse(d1,zero,2.0e-6,2.0e-6);     /* presaturation */ 
  else
    {
     delay(d1_wait);
     rgpulse(d1,zero,2.0e-6,2.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);

  rlpower(tpwr,TODEV);
  rcvroff();
  delay(20.0e-6);

  if (shp_flg[A] == 'y') {

    /* shaped pulse */
    rlpower(tpwrsl,TODEV);
    shaped_pulse(shp_sl,pw_sl,two,2.0e-6,0.0);
    delay(2.0e-6);
    rlpower(tpwr,TODEV); 
    /* shaped pulse */

     delay(2.0e-6);
     rgradient('z',gzlvl1);
     delay(gt1);
     rgradient('z',0.0);
     delay(1000.0e-6);

  }

  rgpulse(pw,t1,1.0e-6,0.0);           /* proton excit. pulse, */

  /* shaped pulse */
  rlpower(tpwrsl,TODEV);
  shaped_pulse(shp_sl,pw_sl,two,2.0e-6,0.0);
  delay(2.0e-6);
  rlpower(tpwr,TODEV); 
  /* shaped pulse */


  txphase(t3);
  dec2phase(t2);

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

  delay(tauxh - gt2 - 4.0e-6 - pw_sl);              /* delay=1/2J(XH)   */

  dec2rgpulse(pwx2,t2,0.0,0.0);

  dec2phase(t4);

  delay(tau1);                  /* delay=t2/2      */ 
  
  rgpulse(2*pw,t3,0.0,0.0); 

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

  dec2rgpulse(pwx2,t4,0.0,0.0);

  /* shaped pulse */
  rlpower(tpwrsl,TODEV);
  shaped_pulse(shp_sl,pw_sl,one,2.0e-6,0.0);
  delay(2.0e-6);
  rlpower(tpwr,TODEV); 
  /* shaped pulse */

  rlpower(dpwr2,DO2DEV);  /* lower decoupler power for decoupling on decouper channel 2 */

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

  delay(tauxh - POWER_DELAY - gt2 - 4.0e-6 - 2.0e-6 - pw_sl);  


  rgpulse(pw,one,2.0e-6,0.0); 

/* acquire data */

status(D);
  setreceiver(t5);
}