Changeset fafb43

Timestamp:

Feb 6, 2009, 9:48:35 AM (17 years ago)

Author:

Frederik Heber <heber@…>

Children:

4aef8a

Parents:

d3d15c

Message:

Purely cosmetical changes

File:

• util/src/NanoCreator.c (modified) (69 diffs)

util/src/NanoCreator.c

@@ -25 +25 @@
  * \return pointer to allocated memory
  */
-void * Malloc (size_t size, const char *msg) 
+void * Malloc (size_t size, const char *msg)
 {
   void *ptr = malloc(size);
@@ -45 +45 @@
  * \return pointer to allocated memory
  */
-void * Calloc (size_t size, double value, const char *msg) 
+void * Calloc (size_t size, double value, const char *msg)
 {
   void *ptr = calloc(size, value);
@@ -76 +76 @@
 
 /** Allocs memory and prints a message on fail.
- * \param **old_ptr pointer to old memory range 
+ * \param **old_ptr pointer to old memory range
  * \param *newsize new size of alloc in bytes
  * \param *msg error msg
  * \return pointer to allocated memory
  */
-void * Realloc (void *old_ptr, size_t newsize, const char *msg) 
+void * Realloc (void *old_ptr, size_t newsize, const char *msg)
 {
   if (old_ptr == NULL) {
@@ -104 +104 @@
  * \return pointer to allocated segment with contents
  */
-char * ReadBuffer (char *filename, int *bufferlength) 
+char * ReadBuffer (char *filename, int *bufferlength)
 {
   if ((filename == NULL) || (bufferlength == NULL)) {
@@ -160 +160 @@
  * \param atomicnumber number of atoms in cell
  */
-void AddAtomicNumber(char *filename, int atomicnumber, double **Vector, double **Recivector) 
+void AddAtomicNumber(char *filename, int atomicnumber, double **Vector, double **Recivector)
 {
   if ((filename == NULL) || (Vector == NULL) || (Recivector == NULL)) {
@@ -175 +175 @@
   double volume = Determinant(Vector);
   time_t now;
-  
+
   now = time((time_t *)NULL);   // Get the system time and put it into 'now' as 'calender time'
   // open for writing and prepend
@@ -204 +204 @@
  * \param *factors pointer to array with length and radius factor
  * \param seed random number seed
- * \param numbercell number of atoms in unit cell, needed as length of \a *randomness 
+ * \param numbercell number of atoms in unit cell, needed as length of \a *randomness
  * \param *randomness for each atom in unit cell a factor between 0..1, giving its probability of appearance
  */
@@ -254 +254 @@
   }
   int i,j;
-  
+
   for (i=0;i<NDIM;i++)
     for (j=0;j<NDIM;j++)
       returnvector[j] += matrixref[i][j] * vectorref[i];
-      
-  return returnvector;  
+
+  return returnvector;
 }
 double *MatrixTrafoInverse(double *vectorref, double **matrixref)
@@ -274 +274 @@
   }
   int i,j;
-  
+
   for (i=0;i<NDIM;i++)
     for (j=0;j<NDIM;j++)
       returnvector[i] += matrixref[i][j] * vectorref[j];
-      
-  return returnvector;  
+
+  return returnvector;
 }
 
@@ -286 +286 @@
  * \param **inverse allocated space for inverse of \a **matrix
  */
-void MatrixInversion(double **matrix, double **inverse) 
+void MatrixInversion(double **matrix, double **inverse)
 {
   if ((matrix == NULL) || (inverse == NULL)) {
@@ -294 +294 @@
   // determine inverse
   double det = Determinant(matrix);
-  inverse[0][0] = (matrix[1][1]*matrix[2][2] - matrix[1][2]*matrix[2][1])/det;  
-  inverse[1][0] = (matrix[0][2]*matrix[2][1] - matrix[0][1]*matrix[2][2])/det;  
-  inverse[2][0] = (matrix[0][1]*matrix[1][2] - matrix[0][2]*matrix[1][1])/det;  
-  inverse[0][1] = (matrix[1][2]*matrix[2][0] - matrix[1][0]*matrix[2][2])/det;  
-  inverse[1][1] = (matrix[0][0]*matrix[2][2] - matrix[0][2]*matrix[2][0])/det;  
-  inverse[2][1] = (matrix[0][2]*matrix[1][0] - matrix[0][0]*matrix[1][2])/det;  
-  inverse[0][2] = (matrix[1][0]*matrix[2][1] - matrix[1][1]*matrix[2][0])/det;  
-  inverse[1][2] = (matrix[0][1]*matrix[2][0] - matrix[0][0]*matrix[2][1])/det;  
+  inverse[0][0] = (matrix[1][1]*matrix[2][2] - matrix[1][2]*matrix[2][1])/det;
+  inverse[1][0] = (matrix[0][2]*matrix[2][1] - matrix[0][1]*matrix[2][2])/det;
+  inverse[2][0] = (matrix[0][1]*matrix[1][2] - matrix[0][2]*matrix[1][1])/det;
+  inverse[0][1] = (matrix[1][2]*matrix[2][0] - matrix[1][0]*matrix[2][2])/det;
+  inverse[1][1] = (matrix[0][0]*matrix[2][2] - matrix[0][2]*matrix[2][0])/det;
+  inverse[2][1] = (matrix[0][2]*matrix[1][0] - matrix[0][0]*matrix[1][2])/det;
+  inverse[0][2] = (matrix[1][0]*matrix[2][1] - matrix[1][1]*matrix[2][0])/det;
+  inverse[1][2] = (matrix[0][1]*matrix[2][0] - matrix[0][0]*matrix[2][1])/det;
   inverse[2][2] = (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0])/det;
-   
+
   // check inverse
-  int flag = 0; 
+  int flag = 0;
   int i,j,k;
-  double tmp;   
+  double tmp;
   fprintf(stdout, "Checking inverse ... ");
   for (i=0;i<NDIM;i++)
@@ -319 +319 @@
         } else {
           flag = (fabs(tmp) > MYEPSILON) ? 1 : 0;
-        }        
+        }
       }
     }
@@ -357 +357 @@
       flip(&matrix[i][j],&matrix[j][i]);
 }
- 
+
 
 /** Computes the determinant of a NDIMxNDIM matrix.
@@ -370 +370 @@
   double det =   matrix[0][0] * (matrix[1][1]*matrix[2][2] - matrix[1][2]*matrix[2][1])
                - matrix[1][1] * (matrix[0][0]*matrix[2][2] - matrix[0][2]*matrix[2][0])
-               + matrix[2][2] * (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0]);  
+               + matrix[2][2] * (matrix[0][0]*matrix[1][1] - matrix[0][1]*matrix[1][0]);
   return det;
 }
@@ -392 +392 @@
   }
   int i;
-  
+
   for (i=0;i<NDIM;i++)
     returnvector[i] = vector1[i] + vector2[i];
-    
+
   return returnvector;
 }
@@ -412 +412 @@
   double betrag;
   int i;
-  
+
   // first vector is untouched
   // second vector
-  betrag = Projection(orthvector[axis[1]], orthvector[axis[0]]); 
+  betrag = Projection(orthvector[axis[1]], orthvector[axis[0]]);
   fprintf(stdout,"%lg\t",betrag);
   for (i=0;i<NDIM;i++)
      orthvector[axis[1]][i] -= orthvector[axis[0]][i] * betrag;
   // third vector
-  betrag = Projection(orthvector[axis[0]], orthvector[axis[2]]); 
+  betrag = Projection(orthvector[axis[0]], orthvector[axis[2]]);
   fprintf(stdout,"%lg\t",betrag);
   for (i=0;i<NDIM;i++)
      orthvector[axis[2]][i] -= orthvector[axis[0]][i] * betrag;
-  betrag = Projection(orthvector[axis[1]], orthvector[axis[2]]); 
+  betrag = Projection(orthvector[axis[1]], orthvector[axis[2]]);
   fprintf(stdout,"%lg\n",betrag);
   for (i=0;i<NDIM;i++)
@@ -433 +433 @@
  * \param *vector1 reference vector
  * \param *vector2 vector which is projected
- * \return projection 
+ * \return projection
  */
 double Projection(double *vector1, double *vector2)
@@ -457 +457 @@
   double scp = 0.;
   int i;
-  
+
   for (i=0;i<NDIM;i++)
     scp += vector1[i] * vector2[i];
-    
+
   return scp;
 }
@@ -468 +468 @@
  * \return norm of \a *vector
  */
-double Norm(double *vector) 
+double Norm(double *vector)
 {
   if (vector == NULL) {
@@ -540 +540 @@
   double diameter[2] = {0.,0.};
   int i, biggest;
-  
+
   for (i=0;i<NDIM;i++) {
     diameter[0] += (matrix[axis[0]][i]*factors[0] - matrix[axis[1]][i]*factors[1]) * (matrix[axis[0]][i]*factors[0] - matrix[axis[1]][i]*factors[1]);
@@ -555 +555 @@
 
   return diameter[biggest];
-} 
+}
 
 /** Determines the center of gravity of atoms in a buffer \a bufptr with given \a number
@@ -576 +576 @@
   char name[255], line[255];
   int i,j;
-  
+
   // Determine center of gravity
   for (i=0;i<number;i++) {
@@ -587 +587 @@
   for (i=0;i<NDIM;i++)
     cog[i] /= -number;
-    
+
   Free(atom, "CenterOfGravity: atom");
   return cog;
@@ -605 +605 @@
   for (i=0; i<NDIM; i++)
     TempVectors[i] = (double *) Malloc(sizeof(double)*NDIM, "Main: TempVectors");
-  
+
   for (i=0; i<NDIM; i++)
     for (j=0; j<NDIM; j++)
@@ -618 +618 @@
   for (i=0; i<NDIM; i++)
     OrthoVector[axis[0]][i] = TempVectors[TempAxis[2]][i]*factor;
-  
+
   TempAxis[1] = axis[1];
   TempAxis[2] = axis[0];
@@ -629 +629 @@
   for (i=0; i<NDIM; i++)
     OrthoVector[axis[1]][i] = TempVectors[TempAxis[2]][i]*factor;
-  
+
   for (i=0; i<NDIM; i++)
     OrthoVector[axis[2]][i] = Vector[axis[2]][i];
@@ -635 +635 @@
   fprintf(stdout, "Orthogonal vectors are: \n");
   for (i=0; i<NDIM; i++) {
-    for (j=0; j<NDIM; j++) 
+    for (j=0; j<NDIM; j++)
       fprintf(stdout, "%lg\t", OrthoVector[axis[i]][j]);
     fprintf(stdout, "\n");
   }
-  
+
   // free memory
   Free(TempAxis, "CreateOrthogonalAxisVectors: *TempAxis");
@@ -664 +664 @@
 // --------------------------------------- M A I N ---------------------------------------------------
 int main(int argc, char **argv) {
-        char *filename = NULL; 
+        char *filename = NULL;
   char *CellFilename = NULL, *SheetFilename = NULL, *TubeFilename = NULL, *TorusFilename = NULL;
   char *SheetFilenameAligned = NULL, *TubeFilenameAligned = NULL;
@@ -678 +678 @@
   int i,j, ggT;
   int length;
-  
-        // Read command line arguments  
+
+        // Read command line arguments
         if (argc <= 2) {
                 fprintf(stdout, "Usage: %s <file> <stage>\n\tWhere <file> specifies a file to start from <stage> or a basename\n\t<stage> is either None, Cell, Sheet, Tube, Torus and specifies where to start the rolling up from.\n\tNote: The .Aligned. files can't be used (rotation is essential).\n", argv[0]);
@@ -688 +688 @@
                 strncpy(stage, argv[2], 5);
     fprintf(stdout, "I will begin at stage %s.\n", stage);
-                
+
                 // build filenames
     char *ptr = strrchr(filename, '.');
     length = strlen(filename);
     if (ptr != NULL) {
-      length = ((int)(ptr-filename) >= length-5) ? (int)(ptr-filename) : length;  
+      length = ((int)(ptr-filename) >= length-5) ? (int)(ptr-filename) : length;
       filename[length] = '\0'; // eventueller
     }
@@ -710 +710 @@
     sprintf(TubeFilenameAligned, "%s.Tube.Aligned.xyz", filename);
         }
-        
+
   // Allocating memory
   Debug ("Allocating memory\n");
@@ -726 +726 @@
     TubevectorInverse[i] = (double *) Malloc(sizeof(double)*NDIM, "Main: TubevectorInverse");
   }
-  
+
   // ======================== STAGE: Cell ==============================
   // The cell is simply created by transforming relative coordinates within the cell
   // into cartesian ones using the unit cell vectors.
-    
+
   double volume;
   int numbercell;
@@ -740 +740 @@
         fprintf(stdout, "You will give the unit cell of the given substance.\nAfterwards, the programme will create a Sheet, a Tube and a Torus, each with their own xyz-file named accordingly.\n\n");
     fprintf(stdout, "Enter 1st primitive vector: ");
-    fscanf(stdin, "%lg %lg %lg", &Vector[0][0], &Vector[0][1], &Vector[0][2]);  
+    fscanf(stdin, "%lg %lg %lg", &Vector[0][0], &Vector[0][1], &Vector[0][2]);
     fprintf(stdout, "Enter 2nd primitive vector: ");
-    fscanf(stdin, "%lg %lg %lg", &Vector[1][0], &Vector[1][1], &Vector[1][2]);  
+    fscanf(stdin, "%lg %lg %lg", &Vector[1][0], &Vector[1][1], &Vector[1][2]);
     fprintf(stdout, "Enter 3rd primitive vector: ");
     fscanf(stdin, "%lg %lg %lg", &Vector[2][0], &Vector[2][1], &Vector[2][2]);
     fprintf(stdout,"Unit vectors are\n");
-    PrintMatrix(stdout, Vector);  
+    PrintMatrix(stdout, Vector);
         } else {
     char *ptr = NULL;
@@ -762 +762 @@
       ptr += strlen(dummy);
       sscanf(ptr, "%lg %lg %lg", &Vector[i][0], &Vector[i][1], &Vector[i][2]);
-      fprintf(stdout, "%5.5lg %5.5lg %5.5lg\n", Vector[i][0], Vector[i][1], Vector[i][2]); 
+      fprintf(stdout, "%5.5lg %5.5lg %5.5lg\n", Vector[i][0], Vector[i][1], Vector[i][2]);
     }
         }
@@ -768 +768 @@
   volume = Determinant(Vector);
   fprintf(stdout,"Volume is %lg\n", volume);
-  MatrixInversion(Vector, Recivector);  
+  MatrixInversion(Vector, Recivector);
   //Transpose(Recivector);   // inverse's got row vectors if normal matrix' got column ones
   fprintf(stdout, "Reciprocal vector is ");
   PrintMatrix(stdout, Recivector);
   fprintf(stdout, "Reciprocal volume is %lg\n", Determinant(Recivector));
-  
+
   fprintf(stdout, "Vector magnitudes: %5.5lg %5.5lg %5.5lg\n", Norm(Vector[0]), Norm(Vector[1]), Norm(Vector[2]));
-   
+
   Debug ("STAGE: Cell\n");
   if (!strncmp(stage, "Non", 3)) {
@@ -793 +793 @@
       fprintf(stdout, "Atom %i: %s %5.5lg %5.5lg %5.5lg\n", i, name, tempvector[0], tempvector[1], tempvector[2]);
       fprintf(stdout, "Probe: %s %5.5lg %5.5lg %5.5lg\n", name,
-        atom[0]*Vector[0][0]+atom[1]*Vector[1][0]+atom[2]*Vector[2][0],      
-        atom[0]*Vector[0][1]+atom[1]*Vector[1][1]+atom[2]*Vector[2][1],      
-        atom[0]*Vector[0][2]+atom[1]*Vector[1][2]+atom[2]*Vector[2][2]    
+        atom[0]*Vector[0][0]+atom[1]*Vector[1][0]+atom[2]*Vector[2][0],
+        atom[0]*Vector[0][1]+atom[1]*Vector[1][1]+atom[2]*Vector[2][1],
+        atom[0]*Vector[0][2]+atom[1]*Vector[1][2]+atom[2]*Vector[2][2]
         );
       fprintf(CellFile, "%s %lg %lg %lg\n", name, tempvector[0], tempvector[1], tempvector[2]);
@@ -805 +805 @@
 
     CellBuffer = ReadBuffer(CellFilename, &length);
-        
+
     sprintf(stage, "Cell");
   } else {
@@ -812 +812 @@
     sscanf(line, "%i", &numbercell);
   }
- 
+
   fprintf(stdout, "\nThere are %i atoms in the cell.\n", numbercell);
-    
+
   // ======================== STAGE: Sheet =============================
   // The sheet is a bit more complex. We read the cell in cartesian coordinates
@@ -838 +838 @@
   FILE *SheetFile = NULL;
   FILE *SheetFileAligned = NULL;
-  
+
   Debug ("STAGE: Sheet\n");
   if (!strncmp(stage, "Cell", 4)) {
     fprintf(stdout, "\nEnter seed unequal 0 if any of the bonds shall have a randomness in their being: ");
     fscanf(stdin, "%d", &seed);
-    if (seed != 0) 
+    if (seed != 0)
       input = 'y';
     randomness = (double *) Malloc(sizeof(double)*numbercell, "Main: at sheet - randomness");
@@ -856 +856 @@
       else { randomness[i] = 1.-percentage; }
     };
-    
+
     fprintf(stdout, "\nSpecify the axis permutation that is going to be perpendicular to the sheet [tubeaxis, torusaxis, noaxis]: ");
     fscanf(stdin, "%d %d %d", &axis[0], &axis[1], &axis[2]);
     fprintf(stdout, "axis: %d %d %d\n", axis[0], axis[1], axis[2]);
-    
+
     // create orthogonal vectors individually for each unit cell vector
     CreateOrthogonalAxisVectors(OrthoVector, Vector, axis);
@@ -866 +866 @@
     fprintf(stdout, "Orthogonal vector axis[1] %lg\n", Projection(Vector[axis[1]], OrthoVector[axis[1]]));
 
-    do {      
+    do {
       fprintf(stdout, "\nNow specify the two natural numbers (m n) defining the chiral angle, \nif the result is crap, try flipping to (m,n): ");
       fscanf(stdin, "%d %d", &chiral[0], &chiral[1]);
@@ -879 +879 @@
         Tubevector[axis[1]][i] = (double)chiral[0] * OrthoVector[axis[0]][i] - (double)chiral[1] * OrthoVector[axis[1]][i];
         //Tubevector[axis[1]][i] = (-chiral[0]-2.*chiral[1])/(double)ggT * Vector[axis[0]][i] + (2.*chiral[0]+chiral[1])/(double)ggT * Vector[axis[1]][i];
-//        fprintf(stderr, "Tubevector[axis[0]][i] = (double)chiral[0] * Vector[axis[0]][i] + (double)chiral[1] * Vector[axis[1]][i]\n = %lg * %lg + %lg * %lg = %lg + %lg = %lg\n\n", 
-//          (double)chiral[0], Vector[axis[0]][i], (double)chiral[1], Vector[axis[1]][i], 
-//          (double)chiral[0] * Vector[axis[0]][i], (double)chiral[1] * Vector[axis[1]][i], 
+//        fprintf(stderr, "Tubevector[axis[0]][i] = (double)chiral[0] * Vector[axis[0]][i] + (double)chiral[1] * Vector[axis[1]][i]\n = %lg * %lg + %lg * %lg = %lg + %lg = %lg\n\n",
+//          (double)chiral[0], Vector[axis[0]][i], (double)chiral[1], Vector[axis[1]][i],
+//          (double)chiral[0] * Vector[axis[0]][i], (double)chiral[1] * Vector[axis[1]][i],
 //          Tubevector[axis[0]][i]);
         Tubevector[axis[2]][i] = Vector[axis[2]][i];
@@ -917 +917 @@
         x1 = gsl_vector_get(u,0)*(double)i;
         x2 = gsl_vector_get(u,1)*(double)i;
-        x3 = 
+        x3 =
         fprintf(stdout, "%d: %d\t%d vs. %lg\t%lg\n",i, ((int)(x1+x1/fabs(x1)*.5)), ((int)(x2+x2/fabs(x2)*.5)), (x1), (x2));
         if (( fabs( ((int)(x1+x1/fabs(x1)*.5)) - (x1) ) < 1e-6) && ( fabs( ((int)(x2+x2/fabs(x2)*.5)) - (x2) ) < 1e-6 )) {
@@ -956 +956 @@
       gsl_vector_free(u);
       gsl_eigen_symmv_free(w);
-      
+
       if (fabs(x1) > 1e-6) {
         fprintf(stderr,"Resulting TubeVectors of axis %d and %d and not orthogonal, aborting.\n", axis[0], axis[1]);
         return(128);
       }
-        
-      
+
+
       angle = Projection(Tubevector[axis[1]], Vector[axis[0]]);
       fprintf(stdout, "Projection Tubevector1 axis[0] %lg %lg\n", angle, 1./angle);
       angle = Projection(Tubevector[axis[1]], Vector[axis[1]]);
       fprintf(stdout, "Projection Tubevector1 axis[1] %lg %lg\n", angle, 1./angle);
-      
+
 /*      fprintf(stdout, "Vector\n");
       PrintMatrix(stdout, Vector);
@@ -1006 +1006 @@
     bufptr += (GetNextline(bufptr, line))*sizeof(char);
     sscanf(line, "%d", &numbersheet);
-    
+
     // retrieve axis permutation
     sprintf(dummy,  "Axis");
@@ -1021 +1021 @@
     ptr += strlen(dummy);
     sscanf(ptr, "%d %d %d", &axis[0], &axis[1], &axis[2]);
-    fprintf(stdout, "%d %d %d\n", axis[0], axis[1], axis[2]); 
-    
+    fprintf(stdout, "%d %d %d\n", axis[0], axis[1], axis[2]);
+
     // retrieve chiral numbers
     sprintf(dummy,  "(n,m)");
@@ -1034 +1034 @@
       fprintf(stderr, "ERROR: Main at stage Sheet - could not find %s in %s\n", dummy, SheetFilename);
       exit(255);
-    }      
+    }
     ptr += strlen(dummy);
     sscanf(ptr, "%d %d", &chiral[0], &chiral[1]);
-    fprintf(stdout, "%d %d\n", chiral[0], chiral[1]); 
+    fprintf(stdout, "%d %d\n", chiral[0], chiral[1]);
     ggT = GCD(2*chiral[1]+chiral[0],2*chiral[0]+chiral[1]);
     fprintf(stdout, "Greatest Common Denominator of (2n+m, 2m+n) is %d\n", ggT);
-    
+
     // retrieve length and radius factors
     sprintf(dummy,  "factors");
@@ -1055 +1055 @@
     ptr += strlen(dummy);
     sscanf(ptr, "%d %d %d", &factors[0], &factors[1], &factors[2]);
-    fprintf(stdout, "%d %d %d\n", factors[0], factors[1], factors[2]); 
+    fprintf(stdout, "%d %d %d\n", factors[0], factors[1], factors[2]);
 
     // create orthogonal vectors individually for each unit cell vector
@@ -1102 +1102 @@
       x1 = gsl_vector_get(u,0)*(double)i;
       x2 = gsl_vector_get(u,1)*(double)i;
-      x3 = 
+      x3 =
       fprintf(stdout, "%d: %d\t%d vs. %lg\t%lg\n",i, ((int)(x1+x1/fabs(x1)*.5)), ((int)(x2+x2/fabs(x2)*.5)), (x1), (x2));
       if (( fabs( ((int)(x1+x1/fabs(x1)*.5)) - (x1) ) < 1e-6) && ( fabs( ((int)(x2+x2/fabs(x2)*.5)) - (x2) ) < 1e-6 )) {
@@ -1140 +1140 @@
     gsl_vector_free(v);
     gsl_vector_free(u);
-    gsl_eigen_symmv_free(w);    
+    gsl_eigen_symmv_free(w);
 
     if (fabs(x1) > 1e-6) {
@@ -1147 +1147 @@
     }
 
-    // retrieve seed ... 
+    // retrieve seed ...
     randomness = (double *) Calloc(sizeof(double)*numbercell, 0., "Main: at sheet - randomness");
     sprintf(dummy,  "seed");
@@ -1159 +1159 @@
       fprintf(stderr, "ERROR: Main at stage Sheet - could not find %s in %s\n", dummy, SheetFilename);
       exit(255);
-    }      
+    }
     ptr += strlen(dummy);
     sscanf(ptr, "%d", &seed);
-    fprintf(stdout, "%d\n", seed); 
+    fprintf(stdout, "%d\n", seed);
 
     // ... and randomness
@@ -1187 +1187 @@
         sscanf(line, "%d %lg", &j, &dummydouble);
         randomness[j] = dummydouble;
-        fprintf(stdout, "%d %g\n", j, randomness[j]);      
+        fprintf(stdout, "%d %g\n", j, randomness[j]);
       }
     }
@@ -1204 +1204 @@
   //angle -= acos(Projection(Vector[axis[0]], Tubevector[axis[0]]));
   //angle = 30./180.*M_PI - acos(Projection(Vector[axis[0]], Tubevector[axis[0]]));
-  angle = -acos(Projection(Tubevector[axis[0]], Vector[axis[0]])); 
+  //angle = acos(Projection(Tubevector[axis[0]], Vector[axis[0]]));
   fprintf(stdout, "The relative alignment rotation angle then is %lg\n", angle/M_PI*180.);
   if (fabs(Tubevector[axis[0]][0]) > MYEPSILON)
-    angle = M_PI/2. - acos(Tubevector[axis[0]][0]/Norm(Tubevector[axis[0]]));
+    angle = -M_PI/2. + acos(Tubevector[axis[0]][0]/Norm(Tubevector[axis[0]]));
   else
     angle = 0.;
@@ -1251 +1251 @@
       sheetnr[j] = sheetnr[j] > tmp ? sheetnr[j] : tmp;
     }
-  } 
+  }
   fprintf(stdout, "Maximum indices to regard: %d %d %d\n", sheetnr[0], sheetnr[1], sheetnr[2]);
   for (i=0;i<NDIM;i++) {
@@ -1260 +1260 @@
     // parse in atoms for quicker processing
     struct Atoms *atombuffer = malloc(sizeof(struct Atoms)*numbercell);
-    bufptr = CellBuffer; 
+    bufptr = CellBuffer;
     bufptr += GetNextline(bufptr, line)*sizeof(char);
     bufptr += GetNextline(bufptr, line)*sizeof(char);
@@ -1301 +1301 @@
         for (nr = 0; nr < numbercell; nr++) {
           percentage = rand()/(RAND_MAX+1.0);
-          //fprintf(stdout, "Lucky number for %d is %lg >? %lg\n", nr, percentage, randomness[nr]); 
+          //fprintf(stdout, "Lucky number for %d is %lg >? %lg\n", nr, percentage, randomness[nr]);
           if (percentage >= randomness[nr]) {
             // Create coordinates at atom site
             coord = VectorAdd(atombuffer[nr].x, offset);
-            //fprintf(stdout, "Atom Nr. %i: ", (numbersheet+1));         
+            //fprintf(stdout, "Atom Nr. %i: ", (numbersheet+1));
             //PrintVector(stdout, coord);
             // project down on major and minor Tubevectors and check for length if out of sheet
@@ -1350 +1350 @@
   //}
   SheetBuffer = ReadBuffer(SheetFilename, &length);
-  
-  
+
+
   // ======================== STAGE: Tube ==============================
   // The tube starts with the rectangular (due to the orthogonalization) sheet
@@ -1365 +1365 @@
   FILE *TubeFile = NULL;
   FILE *TubeFileAligned = NULL;
-  
+
   Debug ("STAGE: Tube\n");
   if (!strncmp(stage, "Sheet", 4)) {
@@ -1385 +1385 @@
     //cog_projected = MatrixTrafoInverse(cog, TubevectorInverse);
     //fprintf(stdout, "\nCenter of Gravity at (%5.5lg\t%5.5lg\t%5.5lg) and projected at (%5.5lg\t%5.5lg\t%5.5lg)\n", cog[0], cog[1], cog[2], cog_projected[0], cog_projected[1], cog_projected[2]);
-    
+
     // restart
     bufptr = SheetBuffer;
@@ -1395 +1395 @@
     for (i=0;i<NDIM;i++)
       halfaxis[i] = factors[0]*Norm(Tubevector[axis[0]])/Norm(Tubevector[i]);
-      
-    double arg, radius;  
+
+    double arg, radius;
     for (i=0;i<numbersheet;i++) {
-      // scan next atom 
+      // scan next atom
       bufptr += GetNextline(bufptr, line);
       sscanf(line, "%s %lg %lg %lg", name, &atom[0], &atom[1], &atom[2]);
@@ -1406 +1406 @@
       //x = VectorAdd(y, cog_projected);
       //free(y);
-      
+
       // roll up (project (x,y,z) on cylindrical coordinates (radius,arg,z))
       arg = 2.*M_PI*x[axis[0]]/(factors[0]) - M_PI;  // is angle
       radius = 1./(2.*M_PI);  // is length of sheet in units of axis vector, divide by pi to get radius (from circumference)
       // fprintf(stdout, "arg: %5.2f (c%2.2f,s%2.2f)\t",$arg, cos($arg), sin($arg));
-      x[axis[0]] = cos(arg)*halfaxis[axis[0]]*radius; //(radius+x[axis[2]]/halfaxis[axis[2]]); // as both vectors are not normalized additional betrag has to be taken into account!
-      x[axis[2]] = sin(arg)*halfaxis[axis[2]]*radius; //(radius+x[axis[2]]/halfaxis[axis[2]]); // due to the back-transformation from eigensystem to cartesian one
+      x[axis[0]] = cos(arg)*halfaxis[axis[0]]*(radius+x[axis[2]]/halfaxis[axis[2]]); // as both vectors are not normalized additional betrag has to be taken into account!
+      x[axis[2]] = sin(arg)*halfaxis[axis[2]]*(radius+x[axis[2]]/halfaxis[axis[2]]); // due to the back-transformation from eigensystem to cartesian one
       //fprintf(stdout, "rotated: (%5.2f,%5.2f,%5.2f)\n",x[0],x[1],x[2]);
       atom_transformed = MatrixTrafo(Tubevector, x);
@@ -1423 +1423 @@
       //fprintf(stdout, "%s\t%5.5lg\t%5.5lg\t%5.5lg\n", name, atom_transformed[0], atom_transformed[1] ,atom_transformed[2]);
 
-      Free(x, "Main: at stage Tube - x"); 
+      Free(x, "Main: at stage Tube - x");
       Free(atom_transformed, "Main: at stage Tube - atom_transformed");
     }
 
-    
+
     fclose(TubeFile);
     fclose(TubeFileAligned);
@@ -1436 +1436 @@
     AddSheetInfo(TubeFilename,axis,chiral, factors, seed, numbercell, randomness);
     fprintf(stdout, "\nThere are %i atoms in the created tube.\n", numbersheet);
-    
+
     strncpy(stage, "Tube", 4);
   } else {
@@ -1454 +1454 @@
 
   FILE *TorusFile;
-  
+
   Debug ("STAGE: Torus\n");
   if (!strncmp(stage, "Tube", 4)) {
@@ -1465 +1465 @@
     bufptr += GetNextline(bufptr, line);  // write numbers to file
     bufptr += GetNextline(bufptr, line);  // write comment to file
-    
+
     //cog = CenterOfGravity(bufptr, numbersheet);
     //cog_projected = MatrixTrafoInverse(cog, TubevectorInverse);
     //fprintf(stdout, "\nCenter of Gravity at (%5.5lg\t%5.5lg\t%5.5lg) and projected at (%5.5lg\t%5.5lg\t%5.5lg)\n", cog[0], cog[1], cog[2], cog_projected[0], cog_projected[1], cog_projected[2]);
-    
+
     // determine half axis as tube vectors not necessarily have same length
     double halfaxis[NDIM];
     for (i=0;i<NDIM;i++)
       halfaxis[i] = Norm(Tubevector[axis[1]])/Norm(Tubevector[i]);
-      
-    double arg, radius;  
+
+    double arg, radius;
     for (i=0;i<numbersheet;i++) {
-      // scan next atom 
+      // scan next atom
       bufptr += GetNextline(bufptr, line);
       sscanf(line, "%s %lg %lg %lg", name, &atom[0], &atom[1], &atom[2]);
@@ -1485 +1485 @@
       //x = VectorAdd(y, cog_projected);
       //free(y);
-      
+
       // roll up (project (x,y,z) on cylindrical coordinates (radius,arg,z))
       arg = 2.*M_PI*x[axis[1]]/(factors[1]) - M_PI;  // is angle
@@ -1497 +1497 @@
       //fprintf(stdout, "%s\t%5.5lg\t%5.5lg\t%5.5lg\n", name, atom_transformed[0], atom_transformed[1] ,atom_transformed[2]);
 
-      Free(x, "Main: at stage Torus - x"); 
+      Free(x, "Main: at stage Torus - x");
       Free(atom_transformed, "Main: at stage Torus - atom_transformed");
     }
-    
+
     fclose(TorusFile);
     //free(cog);
@@ -1536 +1536 @@
   Free(SheetFilenameAligned, "Main: end of stafes - CellFilename");
   Free(TubeFilenameAligned, "Main: end of stafes - CellFilename");
-  
+
   // exit
   exit(0);