/**************************************************************** MSM (Multiplex Spam MQMAS) Au program version: 20050301 ***************************************************************** Authors: C. Fernandez, N. Malicki, and L. Mafra LCS (Laboratoire Catalyse et Spectrochimie) UMR6506 CNRS/ENSICAEN, Universite de Caen-Basse Normandie 6 bd du Marechal Juin, 14050 Caen, France mailto: Christian.fernandez@ensicaen.fr ***************************************************************** Copyright (c)2005 LCS (C.Fernandez, N.Malicki, L.Mafra) This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. ***************************************************************** References: (1) N.Ivchenko et al., J.Magn.Reson. 160 (2003) 52-58 (2) Z.Gan, H.T.Kwak, J.Magn.Reson. 168 (2004) 346-351. (3) N.Malicki, et al., Solid State NMR, 2005 ***************************************************************** Description and Usage: the MSM AU program is designed to process the 3D data generated by the MultiplexSPAMMQMAS pulse program *****************************************************************/ //---------------------------------------------------------------- // Declaration and Initialisation //---------------------------------------------------------------- int *ivector(long nl, long nh); void free_ivector(int *v, long nl, long nh); double *dvector(long nl, long nh); void free_dvector(double *v, long nl, long nh); #define PI 3.1415926535897932384626433832795 //Miscellanous int newexpno, oldexpno; double c, s; double phi; int coherence; double tc, ts; //XWinnmr parameters int parmode; // Number of dimension - 1 int endian; // Bytes order flag int td3; // number of point in a row of the SER file (FID) int td2; // number of row in the SER F2 dimemsion (MQMAS) int td1; // number of row in the SER F1 dimemsion (n phases) // File handling long nbytes; // Number of bytes to read int nbytesread; // Number of lines of TD points actually read FILE *fileSER; // Pointer on the 3D SER file FILE *file2D; // Pointer on the 2D file char nameSER[PATH_LENGTH]; // Name of the original SER file char name2D[PATH_LENGTH]; // Name of the temp SER file TEMP // Pointers on file arrays int *SER; // original SER data (1-3) array int *NewSER; // rebuilded SER data array double *TEMPC0; // cosine part of the temporary ZQ data double *TEMPS0; // sine part of the temporary ZQ data double *TEMPC1; // cosine part of the temporary 1Q data double *TEMPS1; // sine part of the temporary 1Q data //---------------------------------------------------------------- // Initializing Processing //---------------------------------------------------------------- // Get the foreground dataset parameters GETCURDATA; // Check if is a 3D spectra FETCHPARS("PARMODE",&parmode); if ( parmode != 2 ) STOPMSG(" Not a 3D spectrum!"); // Read size of the acquisition (F3) dimension FETCHPARS("TD", &td3); if ( 256*(td3/256) != td3 ) Proc_err (DEF_ERR_OPT, "Please use a TD in the F3 acquisition dimension which is a multiple of 256"); // Read size of the phase (F1) dimension FETCHPAR3S("TD",&td1); // Read size of the MQ (F2) dimension FETCHPAR1S("TD",&td2); // ByteOrda (for compatibilities between O/S) FETCHPARS("BYTORDA",&endian); //---------------------------------------------------------------- // Make a copy of the current dataset and work on this copy // ...Safety belt! //---------------------------------------------------------------- newexpno=expno+4000; GETINT("Enter EXPNO for processing the SER file:", newexpno); // Create a new 2D dataset to save the generated SER file RSER2D("s23",1,newexpno); // Make this new dataset available for further AU statements oldexpno=expno; expno=newexpno; //---------------------------------------------------------------- // Open Files for processing //---------------------------------------------------------------- // Open source file 3D SER and care about opening errors (void)sprintf(nameSER,"%s/data/%s/nmr/%s/%d/ser",disk,user,name,oldexpno); if ((fileSER=fopen(nameSER,"rb"))==NULL) { (void)sprintf(text," I/O Error (Open) \n%s ",nameSER); STOPMSG(text); } // Open destination file 2D SER and care about opening errors (void)sprintf(name2D,"%s/data/%s/nmr/%s/%d/ser",disk,user,name,expno); if ((file2D=fopen(name2D,"wb"))==NULL) { (void)sprintf(text," I/O Error (Open) \n%s ",name2D); STOPMSG(text); } //---------------------------------------------------------------- // Processing of the 3D SER file //---------------------------------------------------------------- //Number of bytes to be read in a (1-3) planes // => number of separated phases (TD1) * size of FID in char (td3) *4 length of INT nbytes=td1*td3*4; // Memory allocation for SER and TEMP arrays SER=ivector(0,td1*td3-1); TEMPC0=dvector(0,td3-1); TEMPS0=dvector(0,td3-1); TEMPC1=dvector(0,td3-1); TEMPS1=dvector(0,td3-1); NewSER=ivector(0,2*td2*td3-1); // Arrays initialization for (i1=0;i1<=td1*td3-1;i1++) { SER[i1]=0; } for (i1=0;i1<=td2*td3*2-1;i1++) { NewSER[i1]=0; } // Set which coherence is selected during the evolution period coherence=3; GETINT("p (coherence selected in the evolution period)?:",coherence); // // Process the TD2 (1-3) planes (MQ evolution) // for (i2=0; i2<=td2-1; i2++) { // Read SER (1-3) plane if ((nbytesread=fread(SER,nbytes, 1,fileSER))<=0) { fclose(file2D); fclose(fileSER); sprintf(text," SER file corrupted!"); STOPMSG(text); } local_swap4((char *)SER,nbytes,endian); // Arrays initialization for (i3=0;i3<=td3-1;i3++) { TEMPC0[i3]=0; // zQ pathway TEMPS0[i3]=0; TEMPC1[i3]=0; // +/-1Q pathways TEMPS1[i3]=0; } // // Process the TD1 phases to separate and recombine the various pathways // This lead to the creation of the Cosine and Sine Part (SRH) // of the resulting SER file // for (i1=0; i1<=td1-1; i1=i1+2){ // Compute the numerical phase phi phi=2.0*PI*(double)(coherence*i1)/(double)td1; c=cos(phi); s=sin(phi); // Contribution from phi3 = 0 for (i3=0;i3<=td3-1;i3++) { tc=c*(double)SER[i3+i1*td3]; ts=s*(double)SER[i3+i1*td3]; TEMPC0[i3]=TEMPC0[i3]+tc; TEMPS0[i3]=TEMPS0[i3]+ts; TEMPC1[i3]=TEMPC1[i3]-ts; TEMPS1[i3]=TEMPS1[i3]+tc; } // Contribution from phi3 = 180° for (i3=0;i3<=td3-1;i3++) { tc=c*(double)SER[i3+(i1+1)*td3]; ts= s*(double)SER[i3+(i1+1)*td3]; TEMPC0[i3]=TEMPC0[i3]-tc; TEMPS0[i3]=TEMPS0[i3]-ts; TEMPC1[i3]=TEMPC1[i3]-ts; TEMPS1[i3]=TEMPS1[i3]+tc; } }// End of loop over TD1 phases for (i3=0; i3 //---------------------------------------------------------------- // useful procedures //---------------------------------------------------------------- int *ivector(long nl, long nh) { // allocate an integer vector with subscript range v[nl..nh] int *v; v=(int *)malloc((size_t) ((nh-nl+2)*sizeof(int))); if (!v) { STOPMSG("allocation failure in vector()"); } return v-nl+1; } void free_ivector(int *v, long nl, long nh) { // deallocate a integer vector allocated with ivector() free((char*) (v+nl-1)); } double *dvector(long nl, long nh) { // allocate a double vector with subscript range v[nl..nh] double *v; v=(double *)malloc((size_t) ((nh-nl+2)*sizeof(double))); if (!v) { STOPMSG("allocation failure in vector()"); } return v-nl+1; } void free_dvector(double *v, long nl, long nh) { // deallocate double vector allocated with dvector() free((char*) (v+nl-1)); } /*****************************************************************/