changes for new data format

Makefile

@@ -5,7 +5,7 @@ MAKEFILE      = Makefile
 VERSION_MAJOR	= 0						# Major Version Number
 VERSION_MINOR	= 2						# Minor Version Number
 ####### Compiler, tools and options
-CUDA			= 1
+CUDA			= 0
 DEFINES       	=  -D_DEBUG
 DEFINES			+= -DhcFloat=float		# floating point precision to be used
 DEFINES			+= -DNUMTHREADS=12		# number of threads to be used in multi-threaded computations

bin/.hide

@@ -0,0 +1 @@
+# file necessary for adding empty directory

doc/Doxyfile

@@ -829,7 +829,7 @@ WARN_LOGFILE           =
 # spaces. See also FILE_PATTERNS and EXTENSION_MAPPING
 # Note: If this tag is empty the current directory is searched.
 
-INPUT                  = ../
+INPUT                  = ./
 
 # This tag can be used to specify the character encoding of the source files
 # that doxygen parses. Internally doxygen uses the UTF-8 encoding. Doxygen uses

engine/hcImage.cpp

@@ -1,8 +1,8 @@
 #include "engine/hcImage.h"
 
-#include "gsl/gsl_rng.h"
-#include "gsl/gsl_randist.h"
-#include "gsl/gsl_cdf.h"
+//#include "gsl/gsl_rng.h"
+//#include "gsl/gsl_randist.h"
+//#include "gsl/gsl_cdf.h"
 
 /*! @param numPixels 			number of pixels in meridional direction
  *  @param northPixelBorder		border of north-most pixel (in sin(latitude)), 1.0 for pixels stretching all the way to the north pole

engine/hcImageFITS.cpp

@@ -1,6 +1,6 @@
 #include "engine/hcImageFITS.h"
 
-#include "fftw3.h"
+//#include "fftw3.h"
 
 void printerror(int status)
 {

engine/hcImageFITS.h

@@ -8,7 +8,7 @@
 #include "src/imageStatistics.h"
 
 #include "fitsio.h"
-#include "gsl/gsl_fft_real.h"
+//#include "gsl/gsl_fft_real.h"
 
 typedef unsigned int uint;
 

main.cpp

@@ -115,7 +115,6 @@ bool parse(int argc, char **argv)
 				else if(command == "help"){			help		 	= true;}
 				else if(command == "rescompr"){		resCompR 		= atoi(optarg);}
 				else if(command == "resmaptheta"){	resMapTheta 	= atoi(optarg);}
-				else if(command == "rss"){			rss				= strtof(optarg, NULL);}
 				else if(command == "resmapphi"){	resMapPhi 		= atoi(optarg);}
 				else if(command == "order"){		order 			= atoi(optarg);}
 				else if(command == "ell")
@@ -125,7 +124,6 @@ bool parse(int argc, char **argv)
 					{
 						printErrMess(__FILE__, __LINE__, "ellipticity must be >= 0.0, you supplied --ell " + string(optarg));
 						abort = true;
-						//return false;
 					}
 					else ellipticity = val;
 				}
@@ -135,7 +133,6 @@ bool parse(int argc, char **argv)
 					{
 						printErrMess(__FILE__, __LINE__, "you supplied --method " + string(optarg) + ", but only 'numeric' and 'shc' are supported");
 						abort = true;
-						//return false;
 					}
 
 					if(!strcmp(optarg, "shc"))			method = METH_SHC;
@@ -149,9 +146,10 @@ bool parse(int argc, char **argv)
 					{
 						printErrMess(__FILE__, __LINE__, "source surface radius must be > 1.0, you supplied --rss " + string(optarg));
 						abort = true;
-						//return false;
 					}
 					else rss = val * r_sol;
+
+					cout << "set rss: " << rss << "\n";
 				}
 				else if(command == "height")
 				{
@@ -160,7 +158,6 @@ bool parse(int argc, char **argv)
 					{
 						printErrMess(__FILE__, __LINE__, "you supplied --height " + string(optarg) + ", but it must hold 1.0 <= height <= r_ss");
 						abort = true;
-						//return false;
 					}
 					else height = val * r_sol;
 				}
@@ -168,7 +165,6 @@ bool parse(int argc, char **argv)
 				{
 					printStdOutMess(__FILE__, __LINE__, "unknown option " + string(long_options[option_index].name) + " with arg " + (optarg ? string(optarg) : "none"));
 					abort = true;
-					//return false;
 				}
 				break;
 
@@ -210,12 +206,11 @@ bool parse(int argc, char **argv)
 		printErrMess(__FILE__, __LINE__, "non-option ARGV-elements: ");
 		while (optind < argc) printErrMess(__FILE__, __LINE__,  string(argv[optind++]) + " ");
 		abort = true;
-		//return false;
 	}
 
 	if(method != METH_NUMERIC && ellipticity != 1.0)
 	{
-		printErrMess(__FILE__, __LINE__, "--ell " + string(ell) + "not supported with computation method " + getStringFromMethodID(method));
+		printErrMess(__FILE__, __LINE__, "--ell " + toStr(ellipticity) + "not supported with computation method " + getStringFromMethodID(method));
 		abort = true;
 	}
 
@@ -281,8 +276,7 @@ int main(int argc, char **argv)
 /*! \mainpage
  *
  *  \brief This PFSS computation suite computes the PFSS solution using photospheric magnetograms from MDI, HMI, GONG, and WSO. It can also create
- *  image-like maps of the magnetic configuration and the expansion factor at arbitrary heights between photsphere and source surface.
- *
+ *  maps of the magnetic configuration and the expansion factor at arbitrary heights between photsphere and source surface.
  *  Commands are given via command line arguments to the PFSS suite. For an example how to use the command-line interface see below.
  *  Several solar observatories are supported and automatically recognized for supplying the photospheric magnetogram.
  *
@@ -296,14 +290,14 @@ int main(int argc, char **argv)
  *
  *  and then build the documentation:
  *
- *		cd PFSS
+ *		cd pfss
  *		make documentation
  *
  *  The documentation will then be produced both as HTML and PDF in the directory PFSS/doc.
  *
  *  # Optimizations (optional)
  *
- *  Several optimizationbs can be adjusted in the Makefile. If a CUDA-capable device is present, setting the variable CUDA to '1' will build the program to employ the CUDA device, which might decrease computation time substantially.
+ *  Several optimizations can be adjusted in the Makefile. If a CUDA-capable device is present, setting the variable CUDA to '1' will build the program to employ the CUDA device, which might decrease computation time substantially.
  *  For this optimization the CUDA runtime environment needs to be installed (please consult NVIDIA's webpage for instructions).
  *
  *  The variable NUMTHREADS limits the number of threads to be utilized for multithreaded execution of the program. For optimal performance it should be the same number as CPU cores in the system.
@@ -319,20 +313,20 @@ int main(int argc, char **argv)
  *
  * 	If successful the PFSS computation suite can be built via
  *
- *		cd PFSS
+ *		cd pfss
  *		make
  *
- *	The binary will be placed in PFSS/bin.
+ *	The binary will be placed in pfss/bin.
  *
  *  # Running the program
  *
- *  If you used make to build the PFSS computation suite, the binary file will be stored in PFSS/bin/.
+ *  If you used make to build the PFSS computation suite, the binary file will be stored in pfss/bin/.
  *
  *  ## Configuration and model output files
  *
  *  Upon execution the binary reads a configuration file, which specifies the data and configuration directories.
- *  The configuration directory contains information about start and stop times of Carrington rotations. The data directory contains all the output from the PFSS computation suite. If you run the binary from the PFSS/bin/ directory without
- *  specifying a configuration file, the default file PFSS/config/config will be used. The default data directoy is then PFSS/data. Please consult this default config file to set your own data directory at a location with enough disk space
+ *  The configuration directory contains information about start and stop times of Carrington rotations. The data directory contains all the output from the PFSS computation suite. If you run the binary from the pfss/bin/ directory without
+ *  specifying a configuration file, the default file pfss/config/config will be used. The default data directoy is then pfss/data. Please consult this default config file to set your own data directory at a location with enough disk space
  *  if the default location is not suitable. Absolute paths in your config file allows the binary to be executed from arbitrary shell locations. Manipulation of the configuration file is only necessary if you have special needs for the location of
  *  the computed solutions.
  *
@@ -354,7 +348,7 @@ int main(int argc, char **argv)
  *
  *	\-\- **map** _filename_
  *
- *	Computes the magnetic configuration at specified height [default: photosphere and source surface]. _filename_ points to the configuration file of the computed solution (ending in \*config.cfg). For additional arguments see below.
+ *	Computes the magnetic configuration at specified height [default: photosphere and source surface]. _filename_ points to the configuration file of the computed solution (ending in \*config.cfg). For additional options see below.
  *
  *	\-\- **batchcompute** _directory_
  *
@@ -369,8 +363,8 @@ int main(int argc, char **argv)
  *	\-\- **rss** _value_		source surface height (multiples of solar radius), _value_ is floating point\n
  *	\-\- **ell** _value_		ellipticity of source surface, _value_ is floating point, [default: 1.0 (spheric)]\n
  *	\-\- **resCompR** _value_	computational grid resolution in radial direction, other directions are determined automatically,
- *								_value_ is unsigned integer
- *	\-\- **method** _value_		solution method to be used for PFSS computation. _value_ is either 'shc' for the classic spherical harmonic coefficient approach or 'numeric' for the finite difference solver.
+ *								_value_ is unsigned integer\n
+ *	\-\- **method** _value_		solution method to be used for PFSS computation. _value_ is either 'shc' for the classic spherical harmonic coefficient approach or 'numeric' for the finite difference solver\n
  *	\-\- **order** _value_		maximum principal order to be used with the SHC approach, _value_ is unsigned integer
  *
  *	**Additional options for map**
@@ -393,11 +387,11 @@ int main(int argc, char **argv)
  *
  *	# Example use cases
  *
- *	All examples are executed from the PFSS/bin-directory, so change there:
+ *	All examples are executed from the pfss/bin-directory, so change there:
  *
- *		cd PFSS/bin
+ *		cd pfss/bin
  *
- *  To perform a PFSS model evaluation with default parameters for the synoptic photospheric magnetogram found at PFSS/data/input/synop_Ml_0.2066.fits:
+ *  To perform a PFSS model evaluation with default parameters for the synoptic photospheric magnetogram found at pfss/data/input/synop_Ml_0.2066.fits (not included in the repository, you need to download the magnetogram and place it at that location):
  *
  *  	./pfss --compute ../data/input/synop_Ml_0.2066.fits
  *
@@ -407,9 +401,10 @@ int main(int argc, char **argv)
  *
  *	To generate a magnetic mappings of the same solution at height r=2.4 r_sol with a resolution of 130 x 200 pixels:
  *
- *		./pfss --map ../data/2066/2066_MDI_Kiel_PFSS2.50_GridSph_35x87x175_config.cfg --height 2.4 --resMapTheta 130 --resMapPhi 200
+ *		./pfss --map ../data/2066/2066_MDI_Kiel_PFSS2.50_GridSph_35x87x175_config.cfg --height 2.4 \
+ *		--resMapTheta 130 --resMapPhi 200
  *
- *	To evaluate the PFSS model for all photospheric magnetic maps found in directory PFSS/data/batchinput with a radial grid resolution of 40 grid points and a source surface radius of r=3.1 r_sol:
+ *	To evaluate the PFSS model for all photospheric magnetic maps found in directory pfss/data/batchinput with a radial grid resolution of 40 grid points and a source surface radius of r=3.1 r_sol:
  *
  *		./pfss --batchcompute ../data/input/ --resCompR 40 --rss 3.1
  */

src/magMapping.cpp

@@ -333,8 +333,8 @@ bool MagMapping::exportBinary(const string &filename)
     bool retval 		= true;
     uint sizeofFloat 	= sizeof(hcFloat);
     std::ofstream file;
-    file.open(filename, std::ios::out | std::ios::binary);
-
+    file.open(filename, std::ios::out);//std::ios::binary
+    
     file.write(reinterpret_cast<char*>(&sizeofFloat),	sizeof(uint));
     retval &= info.exportBinary(file);
     file.write(reinterpret_cast<char*>(&sinLatFormat),	sizeof(bool));
@@ -343,12 +343,26 @@ bool MagMapping::exportBinary(const string &filename)
     file.write(reinterpret_cast<char*>(&numTheta),		sizeof(uint));
     file.write(reinterpret_cast<char*>(&numPhi),		sizeof(uint));
     file.write(reinterpret_cast<char*>(&height),		sizeofFloat);
-
+    
+    /*    file<<sizeofFloat<<" ";
+    retval &= info.exportBinary(file);
+    file<<sinLatFormat<<" ";
+    file<<maxSinLat<<" ";
+    file<<compCoords<<" ";
+    file<<numTheta<<" ";
+    file<<numPhi<<" ";
+    file<<height<<" ";
+    
+    file<<std::endl;//remove this again*/
     for(uint i=0; i<numTheta*numPhi; ++i)
     {
-    	file.write(reinterpret_cast<char*>(&coords[i](0)),sizeofFloat);
-        file.write(reinterpret_cast<char*>(&coords[i](1)),sizeofFloat);
-        file.write(reinterpret_cast<char*>(&coords[i](2)),sizeofFloat);
+      file.write(reinterpret_cast<char*>(&coords[i](0)),sizeofFloat);
+      file.write(reinterpret_cast<char*>(&coords[i](1)),sizeofFloat);
+      file.write(reinterpret_cast<char*>(&coords[i](2)),sizeofFloat);
+	    /*file<<coords[i](0)<<" ";
+      file<<coords[i](1)<<" ";
+      file<<coords[i](2);
+      file<<std::endl;//remove this again*/
     }
 
     for(uint j=0; j<numTheta; ++j)
@@ -366,7 +380,7 @@ bool MagMapping::importBinary(const string &filename)
 
     PFSSsolutionInfo inf;
     ifstream stream;
-	stream.open(filename, std::ios::out | std::ios::binary);
+    stream.open(filename, std::ios::out); //std::ios::binary
     bool retval 		= true;
     uint sizeoffloat	= 0;
 	stream.read(reinterpret_cast<char*>(&sizeoffloat),	sizeof(uint));
@@ -379,7 +393,7 @@ bool MagMapping::importBinary(const string &filename)
     stream.read(reinterpret_cast<char*>(&numTheta),		sizeof(uint));
     stream.read(reinterpret_cast<char*>(&numPhi),		sizeof(uint));
     stream.read(reinterpret_cast<char*>(tempFloat),		sizeoffloat); height 	= sizeoffloat==4 ? *(reinterpret_cast<float*>(tempFloat)) : *(reinterpret_cast<double*>(tempFloat));
-
+    //    stream.readline();
     init(inf, sinLatFormat, compCoords, maxSinLat, numTheta, numPhi, height);
 
     for(uint i=0; i<numTheta*numPhi; ++i)
@@ -387,6 +401,7 @@ bool MagMapping::importBinary(const string &filename)
         stream.read(reinterpret_cast<char*>(tempFloat),	sizeoffloat); coords[i](0) = sizeoffloat==4 ? *(reinterpret_cast<float*>(tempFloat)) : *(reinterpret_cast<double*>(tempFloat));
         stream.read(reinterpret_cast<char*>(tempFloat),	sizeoffloat); coords[i](1) = sizeoffloat==4 ? *(reinterpret_cast<float*>(tempFloat)) : *(reinterpret_cast<double*>(tempFloat));
         stream.read(reinterpret_cast<char*>(tempFloat),	sizeoffloat); coords[i](2) = sizeoffloat==4 ? *(reinterpret_cast<float*>(tempFloat)) : *(reinterpret_cast<double*>(tempFloat));
+	//	stream.readline();
     }
 
     for(uint j=0; j<numTheta; ++j)
@@ -1010,7 +1025,7 @@ bool MagMapping::exportASCII(string filename, hcFloat *heights, uint numHeights,
 		printErrMess(__FILE__, __LINE__, "requested file '" + filename + "' cannot be created");
 		return false;
 	}
-
+	std::cout<<"asci map "<<filename.c_str()<<std::endl;
     FILE *output = fopen(filename.c_str(), "w");
 
     char outputstring[100000];

src/pfssSolution.cpp

@@ -387,7 +387,7 @@ bool PFSSsolution::mapHeightLevel(hcFloat height, hcFloat maxSinLat, uint numThe
 	bool sinLatGrid			= numTheta==0 ? solver.grid->sinLatGrid	: sinLatFormat;
 
 	string filenameMap		= getFilename_magMappingBin(info, height, sinLatGrid, compCoords, numTheta_m, numPhi_m);
-
+	string filenameMapAscii         = getFilename_magMappingASCII(info, height, sinLatGrid, compCoords, numTheta_m, numPhi_m);
 	MagMapping map;
 	if(doesFileExist(filenameMap))
 	{
@@ -403,6 +403,7 @@ bool PFSSsolution::mapHeightLevel(hcFloat height, hcFloat maxSinLat, uint numThe
 		printStdOutMess(__FILE__, __LINE__, "map at height " + toStr(height) + " m does not exist, start computation.");
 		map.createAtHeight(info, solver, height, numTheta_m, numPhi_m, maxSinLat_m, sinLatGrid, compCoords);
 		map.exportBinary(filenameMap);
+		map.exportASCII(filenameMapAscii, NULL, 0, solver.grid->lowerR, solver.grid->upperR);
 	}
 
 	map.exportImagePolarity();

src/pfssSolution_batches.cpp

@@ -14,7 +14,7 @@ using namespace filesystem;
 extern string dirData;
 
 /*! @param	inDir					direcory to be worked upon
- * 	@param	r_ss					heliocenric position of source surface
+ * 	@param	r_ss					heliocenric position of source surface (m)
  *	@param	resCompR				radial resolution of computational grid
  *	@param	order					maximum principal order of SHC approach
  */
@@ -57,7 +57,7 @@ bool PFSSsolution::batchKielSHC(const string &inDir, hcFloat r_ss, uint resCompR
  *  computational work, for the GUI version please see HelioMagfield::batchRotationComputing
  *
  *	@param	inDir					direcory to be worked upon
- *	@param	r_ss					heliocenric position of source surface
+ *	@param	r_ss					heliocenric position of source surface (m)
  *	@param	resCompR				radial resolution of computational grid
  *	@param	ellipticity				ellipticity of source surface (default: 1.0 - spherical)
  *