/*    zusammen.c      Programm haengt blockweise gemessene Daten aus HASYLAB
                      aneinander zu einer DAT Datei fuer SIMREF   

Uebersetzen mit: cc -lm zusammen.c  (sonst findet es die Mathe-Funktionen nicht)

Eingabe: * Dateinamen der Eingabedatei im HASYLAB Datenformat.
         Sie enthaelt Abfolge von Bloecken in der richtigen 
         Reihenfolge zu groesseren 2 theta Winkeln hin.         
         Wenn die Bloecke keine gleiche 2-theta Schrittweite besitzen,
         koennen sie dennoch eingeordnet werden!
         * Eine Liste mit Korrekturfaktoren fuer die Intensitaten je eines
         Blocks (diese Korrekturfaktoren sind normalerweise aus der
         Messung von Eichreflexen zwischen den Bloecken gewonnen).
         Format der Liste: Blocknr., ungefaehres 2theta des Reflexes, Faktor.
Ausgabe: * Eine Datei mit nur einem Block, mit den aneinandergehaengten
         Messwerten (xxxxx.zus) im DAT-Format: 
               2theta, Messwert, sigma(Messw), Untergrundwert.
         Beliebige 2theta Werte in aufsteigender Reihenfolge koennen 
         vorkommen, eine feste Schrittweite muss nicht eingehalten werden.

         * Eine Datei mit den EXCLUDED REGIONS (xxx.log), dass heisst,
         den Winkelbereichen, die garnicht gemessen wurden.
  
Erstellt durch K. Burger, Okt. 1994 (zusammen.c), modifiziert fur DAT
statt DUO Dateien Sept. 1996.

*/

/* >>>>>>>>>>>>>>   I N C L U D E -- D A T E I E N   <<<<<<<<<<<<<<<<<< */

#include "/usr/include/stdio.h"
#include "/usr/include/stdlib.h"
#include "/usr/include/math.h"
#include "/usr/include/string.h"
#define MAXWERTE 1000
#define MAXGESWERTE 20000
#define MAXBLOCK 130

int redeviel=0;     /* nur zum Testen */
double banfang,bende,bsw, bmonitor;   /* extern */

/* --------------------------------------------------------------- 
 Liest einen Datenblock im B2 HASYLAB Format 
 Titelzeilen sind schon weg. zuerst Kopfzeile des Blocks, am
 Ende dann 0000.0 <RET> 
--------------------------------------------------------------- */

long int blockles(FILE *in, float *px)    /* px zeigt auf x[0] */

{ int ir,schluss;
  long int n,n1,i;
  char buffer[250];
  char *form1 = "%f %lf %lf %lf %i",
       *form2 = "%f %f %f %f";
  float temp, wert, mon, zeit, winkel;
  int blocknr;

  banfang=bsw=bende=0; 

/* erste Zeile enthaelt "b " und 5 Worte   */

  fgets(buffer,250,in);
  if ( feof(in) ) return 0 ;
  if ( strncmp(buffer,"b temp",6)!=0 ) return -1 ;

/* zweite Zeile enthaelt Blockdaten  */

  fgets(buffer,250,in);
  ir = sscanf(buffer,form1,&temp,&banfang,&bsw,
              &bende,&blocknr); 
  if ( (banfang >= bende) || (banfang<=0) || (bsw==0) || (ir<5))
    {fprintf(stderr,"Blockzeile seltsam! \n >>%s<< ",buffer);
     fprintf(stderr,"Von %f mit %f bis %f \n",banfang,bsw,bende);
     return -1;}
  n = ((bende-banfang)/bsw +.00001) + 1;
  if (n>MAXWERTE) {printf("Zuviele Messpunkte in Block!!!\n");
                   return -1; }
  if (redeviel)
    fprintf(stderr,"Von %f mit %.4f bis %f, n=%d \n",banfang,bsw,bende,n);

/* n Messwerte einlesen */

  n=0; schluss=0;
  while (!schluss && !feof(in))
  {
   fgets(buffer,250,in);
   if (0==strncmp(buffer,"      0.0000",12) )
      schluss=1; /* Block zuende, eventuell verfr"uht  */
   else
     {
     ir = sscanf(buffer,form2,&wert,&mon,&zeit,&winkel);
     if (redeviel)
        printf("%.0f %.1f %.3f \n ",wert,mon,winkel);
     if (  feof(in) ||  ir<4 )
       { printf("\n Fehler beim einlesen in <blockles> in Zeile ");
	 printf("%d\n",i);
	 printf("%s",buffer);
	 schluss=1; }
     else
       { *++px = wert;
	 n++; }
     } /* else */
  }
  bende=banfang+(n-1)*bsw; bmonitor=mon;
  return n;
}
/* --------------------------------------------------------------- */
/*  Schreiben eines Blockes in ein DUO-Datei                       */
/* --------------------------------------------------------------- */

void blockschreib(FILE *out, float y[], float start, float sw, float ziel,         		  long int n)
{ long int i;
  if (n>0)
  {
   /* Schreiben der Blockbefehlszeile fuer das DUO-File */
   fprintf(out,"B\n       0.000  %10.5f  %9.5f   %10.5f    %8.4f  %2d   %4d\n"	          ,start,sw,ziel,300.0,1,n);
   /* Schreiben der Daten */
   for (i=1; i<=n; i++)
    { fprintf(out,"%7.0f ",y[i]);
      if ((i % 8)  <= 0)
        fprintf(out,"\n");
    }
   fprintf(out,"\n");
  }
}

/* --------------------------------------------------------------- */
/*  Hauptprozedur 
    Erzeugt aus einem HASYLAB-Datei mit mehreren Bloecken 
    eine DAT-Datei, wobei die Luecken nicht aufgefuellt werden.  
    Ausserdem wird eine Protokolldatei geschrieben (.log).
    Gleichzeitig wird eine Datei mit Korrekturfaktoren eingelesen, mit
    denen die einzelnen Messwerte der Bloecke multipliziert werden: pro
    Block ein Korr-Faktor.*/ 
/* --------------------------------------------------------------- */

void konvertiere(char *quelle,char *z1datei,char *z2datei)
{ FILE *in, *out1, *out2, *korr;
  fpos_t filepos;
  char buffer[250],korrfile[50];
  double anfang, sw, ende, a, b, fdiff;
  float bwerte[MAXBLOCK][MAXWERTE], banf[MAXBLOCK], bschluss[MAXBLOCK],
        bschritt[MAXBLOCK], bmon[MAXBLOCK], x[MAXWERTE], y;
  float korrtheta[MAXBLOCK], korrfaktor[MAXBLOCK];
  int bmarkiert[MAXBLOCK], korrnr[MAXBLOCK];
  long int bnn[MAXBLOCK];
  long int nges, n=0, nn, i, j, k, iret, idiff;
  int ib, nblock, ikorr, lkorr, nkorr;

  /*  initialisiere markierungsfeld fuer zu ignorierende Bloecke */  
  for (i=1; i<=MAXBLOCK; i++) bmarkiert[i]=0;

  printf("\nGebe Namen der Datei mit den Korrfaktoren: ");
  scanf("%s",korrfile);
  printf("%s\n",korrfile);
  
  if ( ((in=fopen(quelle,"rt")) == NULL) || feof(in) )
   {  fprintf(stderr, "Cannot open input file %s.\n",quelle);
      exit(1); }
  if ( (out1=fopen(z1datei,"wt")) == NULL )
   {  fprintf(stderr, "Fehler beim Oeffnen der Datei %s.\n",z1datei); 
      exit(1); 
   }
  if ( (out2=fopen(z2datei,"wt")) == NULL )
   {  fprintf(stderr, "Fehler beim Oeffnen der Datei %s.\n",z2datei); 
      exit(1);
   }
  if ( (korr=fopen(korrfile,"rt")) == NULL )
   {  fprintf(stderr, "Fehler beim Oeffnen der Datei %s.\n",korrfile); 
      exit(1);
   }

/* Lesen der Korrfaktoren */

   fgets(buffer,250,korr);
   fgets(buffer,250,korr);
   fgets(buffer,250,korr);
   i=0; 
   while (!feof(korr))
     {fgets(buffer,250,korr);
      if (!feof(korr))
      { i++;
        sscanf(buffer,"%d %f %f",&korrnr[i], &korrtheta[i], &korrfaktor[i]);
/*        printf("%3d %.2f %.4f\n",korrnr[i], korrtheta[i], korrfaktor[i]);
*/
      }
     }
   fclose(korr); nkorr=i;
   printf("%3d Korrekturfaktoren eingelesen.\n", nkorr);

/* Umkopieren der Titelzeilen */

  fprintf(out1,"Datei %s konvertiert mit z.c\n",quelle);
  i=0;
  while (!feof(in) && (i++<40))
   {  fgets(buffer,250,in);
      if (0==strncmp(buffer,"titelende",9))
         i=42;
      else 
       { fprintf(out1,"%s",buffer);fprintf(out2,"%s",buffer);} 
   }
  fprintf(out2,"Datei %s konvertiert mit zusammen.c\n",quelle);
  

/* ================== Hauptteil =====================================   */

/* -------- Einlesen aller Daten-Bloecke  --------------------------- */    

  /* anfang, ende,nges      gehoeren zum gesamten Datenblock,
     banfang,bsw,bende      gehoeren zum "neuen" Block      */

  nges=0; anfang=ende=0;
  ib=0; lkorr=0; 
  while (!feof(in))
     {
     /* Einlesen des Datenblocks, Ergebnis skaliert in bwerte[]  */
     iret = blockles(in,x); 
     if (iret<0)
       { printf(
          "\n Schlimmer Fehler beim Einlesen der Daten in blockles \n"); 
         if (feof(in))
           printf("End of File gefunden\n");
         break;
       }
     if ((iret<5)) continue; 

     /* nun hat der Block also >0 Werte */
     ib++;

     if (nges > MAXGESWERTE)
       {printf("Maximale Messpunktezahl ueberschritten -- Abbruch.");
        exit(1);}
     if ( (ib >= MAXBLOCK) ) 
       {printf("Zuviele Bloecke -- Abbruch. ");
        exit(1);}

     /* Korrfaktor suchen */
     /* Falls ein leerer Block in der Datendatei garnicht in die hiesige
        Blockliste aufgenommen wurde, steht der richtige Korr-faktor etwas
        spaeter in der Liste und wird anhand des 2theta Wertes entdeckt. 
        lkorr enthaelt die Nummer des letzten verwendeten Korrfaktors */
     y=-1000;
     for (i=lkorr+1; (i<=lkorr+4) && (i<=nkorr); i++) 
       { if ( (korrtheta[i]>=banfang) && (korrtheta[i]<=bende) ) 
          { y=korrfaktor[i]; lkorr=i;  break;}
       };
     if (y<-999) { printf(
                  "Fehler: fuer Block %3d keinen Korrfaktor gefunden\n", ib); 
                  exit(1); }
     if (y < -0.00001) 
       {printf(">> Block wird ignoriert\n"); ib--; continue;}
     else 
        printf("Korrfaktor Nr %3d (%.4f) fuer Block %3d\n",i,y,ib);

     /* block in Liste merken */
     n=iret; nges+=n;
     bnn[ib]=n; bschritt[ib]=bsw; bmon[ib]=bmonitor;
     banf[ib]=banfang; bschluss[ib]=bende;
     for (i=1; i<=n; i++)
       bwerte[ib][i]=x[i]*y;

     printf("Block Nr %d von %.3f bis %.3f mit sw=%.3f\n",ib,banfang,bende,bsw);
     fprintf(out2,
            "Block Nr %d von %.3f bis %.3f mit sw=%.3f\n",ib,banfang,bende,bsw);
     if ((banfang<anfang)||(ib==1)) anfang=banfang;
     if ((bende>ende)||(ib==1)) ende=bende;

    } /* Ende der While-Schleife ueber Eingabe-Datei */
/* ----------------------------------------------------- */
     nblock=ib;
     fclose(in);

/*   nun enthalten 
     banf[], bschluss[]: Blockanfang 2theta, Blockende,
     bwerte[ib][ix]:      Messwerte
     bschritt[]:         Schrittweite,
     bnn[]:              Anzahl der Messwerte in diesem Block
     bmon[] :            Monitorcounts/Schritt.
     anfang, ende:       kleinster und groesster Winkel des Blocks
     nges:		 gesamte Zahl an Messpunkten in allen Bloecken.
     nblock:             gesamte Zahl von Messbloecken.

     Nun muessen doppelte Bloecke gemittelt werden, mit entsprechend
     kleineren sigmas.

   -------- suche zusammenhaengenden Bereich in 2 theta ---------------       

*/
     printf("\nGesamter Datensatz: Anfang=%.3f  Ende=%.3f   Messpunkte=%6d\n",
            anfang,ende,nges);
     fprintf(out2,
            "Gesamter Datensatz: Anfang=%.3f  Ende=%.3f   Messpunkte=%6d\n",
            anfang,ende,nges);
     if (nges==0) {printf("Keine Daten gefunden.\n"); exit(1);}

     i=0; ib=2; 
     nn=bnn[1]; banfang=banf[1]; bende=bschluss[1];

     while ((i==0) && (ib <= nblock))
     {if (banf[ib]>=bende) 
        { /* ok:  hoeherer Winkel */
         nn+=bnn[ib]; 
         bende=bschluss[ib]; 
         ib++;  }
      else {i=-1;}  /* nicht ok: hier abbrechen */ 
      } 
     ib--;
     printf("Gefundener aufeinanderfolgender 2theta-Bereich:\n");
     printf("Anfang=%.3f  Ende=%.3f   Messpunkte=%6d \n",
            banfang,bende,nn);
     printf("Letzter Block: %4d \n",ib);
     fprintf(out2,"Gefundener aufeinanderfolgender 2theta-Bereich:\n");
     fprintf(out2,"Anfang=%.3f  Ende=%.3f   Messpunkte=%6d \n",
            banfang,bende,nn);
     fprintf(out2,"Letzter Block: %4d \n",ib);

/* -------- merge mehrfache Bloecke ------------------------------------ */    
     if (ib < nblock ) 
     { printf(
       "Die folgenden Bloecke werden auf Doppelte im Messbereich abgesucht\n"); 
       fprintf(out2,
       "Die folgenden Bloecke werden auf Doppelte im Messbereich abgesucht\n"); 
       for (i=(ib+1); i<=nblock; i++)   /* Schleife ueber restliche Bloecke */
       { printf(" %3d ",i);
         for (j=1; j<=ib; j++)          /* Schleife ueber erste Bloecke */
         { if ( (banf[i]==banf[j])&&(bschluss[i]==bschluss[j])
             &&(bschritt[i]==bschritt[j])  )
           { /* dieser Block ist gleich! Also dazuaddieren bzw. mergen */
             printf("Addiere Block %3d zu %3d\n", i, j);
             fprintf(out2,"Addiere Block %3d zu %3d\n", i, j);
             if ( (bmon[i]<=0) || (bmon[j]<=0) )
               {printf("Monitorcounts=0!!! -- Abbruch."); exit(1);}
             for (k=1; k<=bnn[j]; k++)
                bwerte[j][k] += (bwerte[i][k]*bmon[j]/bmon[i]) ;         
             ++bmarkiert[j];  /* flag setzen, dass hier dazuaddiert wurde */ 
             break; /* aus innnerer for schleife */
            } /* end if gleich */
         } /* end for j */
       }  /* end for i */
     }   /* end if */

/* -------- 3. Schreibe Bloecke fortlaufend in DAT Datei ---------------- */    

   fprintf(out1,"DAT  10 %.4f %.4f %.4f   290  1 %5d\n",
                banfang,0.002,bende,nn);
   fprintf(out2,"EXCLUDED REGIONS\n");
   printf("\n\nSchreiben in Datei: \n");
  
   for (i=1; i<=ib; i++)
   { /* Schleife ueber alle zusammenhaengenden Bloecke */
     printf(" %3d ",i);
     if (bmarkiert[i]==0) 
       {for (j=1; j<=bnn[i]; j++)
          fprintf(out1,"%.3f %.1f %.3f %.1f \n", banf[i]+bschritt[i]*(j-1),
                  bwerte[i][j], sqrt(bwerte[i][j]), 0); 
       }
     else
       {for (j=1; j<=bnn[i]; j++)
         {n = bmarkiert[i]+1; y = bwerte[i][j]/n;
          fprintf(out1,"%.3f %.1f %.3f %.1f \n", banf[i]+bschritt[i]*(j-1),
                  y, sqrt(y/n), 0);
         }
       } /* end else */  

      /* excluded region schreiben */ 
      if (i>1) fprintf(out2,"EXC> %8.3f %8.3f \n",
          bschluss[i-1]+bschritt[i-1],banf[i]-bschritt[i]);

     } /* ende schleife alle bloecke */

   fclose(out1); fclose(out2);
} /* Ende konvertiere */


/* -------------------------------------------------------------------- */
/* >>>>>>>>>>>>>>>>>>   H A U P T P R O G R A M M   <<<<<<<<<<<<<<<<<<< */
/* -------------------------------------------------------------------- */

void main(int argc, char *argv[])
{  
   char quelle[22],z1datei[22],z2datei[22];
   int i;

   if (argc != 2)
     { printf("Bitte einen Dateinamen beim Aufruf angeben.\n"); 
       exit(0); }
   else
   {
     strcpy(quelle,argv[1]);
     printf("Quelle %s ---> ",quelle);
     /* nun Namen bilden mit neuer Endung .duo und .log */
     i=0;
     while ((quelle[i] != '.') && (quelle[i] != '\0'))
     { z1datei[i] = z2datei[i] = quelle[i];
       i++;
     }
     z1datei[i] = z2datei[i] = '.';
     z1datei[i+1] = 'z'; 
     z1datei[i+2] = 'u'; 
     z1datei[i+3] = 's'; 
     z2datei[i+1] = 'l';
     z2datei[i+2] = 'o';
     z2datei[i+3] = 'g';
     z1datei[i+4] = z2datei[i+4] = '\0';

     printf("Ziel: %s \n",z1datei);

     konvertiere(quelle,z1datei,z2datei);  
     printf("\n FERTIG. \n");
   }
}