00001 #include "parametermap.h"
00002 #include "cheezyparser.h"
00003 #include "oscar.h"
00004 #include "oscar_source_generator_1d.h"
00005 #include "oscar_correlation_generator_1d.h"
00006 #include "oscar_source_generator_3dcart.h"
00007 #include "oscar_source_generator_3dsphr.h"
00008 #include "oscar_correlation_generator_3dcart.h"
00009 #include "message.h"
00010 #include "basisfunc_imager1d.h"
00011 #include "expander.h"
00012 #include "uncoupled_imager3d.h"
00013 #include "corr3d_ylm.h"
00014 #include "tnt_array1d.h"
00015 #include <iostream>
00016 #include <iomanip>
00017
00018
00019
00020
00021
00022
00023 CSourceFtn1dHisto getSource1d( vector<COSCARLine> particleList, parameterMap m );
00024 CSourceFtn1dHisto getSource1dOldWay( vector<COSCARLine> particleList, parameterMap m );
00025 CSourceFtn3dHisto getSource3dCart( vector<COSCARLine> particleList, parameterMap m );
00026 CSourceFtn3dSphr< CSourceFtn1dHisto >
00027 getSource3dSphr( vector<COSCARLine> particleList, parameterMap m );
00028
00029
00030 CCorrFtn1dHisto getCorrelation1d( vector<COSCARLine> particleList, parameterMap m );
00031 CCorrFtn1dHisto getCorrelation1dOldWay( const CSourceFtn1dHisto& souin, parameterMap m );
00032 CCorrFtn3dHisto getCorrelation3dCart( vector<COSCARLine> particleList, parameterMap m );
00033 CCorrFtn3dSphr getCorrelation3dSphr( vector<COSCARLine> particleList, parameterMap m );
00034
00035
00036 void plot3dSlices( CObject3d& obj, string xFileName, string yFileName, string zFileName );
00037 void plot3dSlices( CHistogram3d& obj, string xFileName, string yFileName, string zFileName );
00038 void plot1d( CObject1d& obj, string outFileName );
00039 void plot1d( CHistogram1d& obj, string outFileName );
00040
00041
00042 void make1dSourceAndCorrelation( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots );
00043 void make1dSourceAndCorrelationOldWay( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots );
00044 void make3dCartSourceAndCorrelation( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots );
00045 void make3dCartExpandedSourceAndCorrelation( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots );
00046 void make3dSphrSourceAndCorrelation( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots );
00047 void makeCRABCorrelations( vector<COSCARLine> particleList, parameterMap& inMap, bool make_plots );
00048 void makeSourceMakerSources( vector<COSCARLine> particleList, parameterMap& inMap, bool make_plots );
00049
00050
00051
00052
00053
00054
00056 void getHelp(void){
00057 cout<<"\nUsage: crab <mode> [inputFile.dat]"<<endl;
00058 cout<<" Valid modes are: -1d : Generate 1d source and correlation"<<endl;
00059 cout<<" -1doldway : Generate 1d source and correlation, using old code for source"<<endl;
00060 cout<<" -3dcart : Generate 3d source and correlation directly in Cartesian"<<endl;
00061 cout<<" -3dcart_ex : Generate 3d source and correlation directly in Cartesian, then expand"<<endl;
00062 cout<<" -3dsphr : Generate 3d source and correlation directly in Ylm's"<<endl;
00063 cout<<" -crab : Generate 1d & 3d correlations and output them in the style of Scott Pratt's CRAB code"<<endl;
00064 cout<<" -sourcemaker : Generate 1d source & 3d source slices and output them in the style of Dave Brown's SourceMaker code"<<endl;
00065 cout<<" -skip_correlation : Skip generation of the correlations if accompanying modes would generate them"<<endl;
00066 cout<<" -make_plots : Make plots for all the sources and correlations generated"<<endl;
00067 cout<<" -h, -help, --help : Prints this message, then quits"<<endl;
00068 exit(0);
00069 }
00070
00071
00072 int main(int argc, char* argv[]){
00073
00074 cout << "*** sample 1d implementation of the CoRrelation AfterBurner (CRAB) using CorAL ***"<<endl;
00075 cout << endl;
00076
00077 MESSAGE << CMessage::warning;
00078
00079
00080 bool skip_correlation=false;
00081 bool make_plots=false;
00082
00083
00084 if (argc==1) getHelp();
00085 string paramFile("");
00086 vector<string> modeList;
00087 for (int iarg = 1; iarg<argc; ++iarg){
00088 string sarg(argv[iarg]);
00089 if (sarg=="-help") getHelp();
00090 if (sarg=="--help") getHelp();
00091 if (sarg=="-h") getHelp();
00092 if (sarg=="-skip_correlation") {
00093 cout << "Will skip correlation function generation\n"<<endl;
00094 skip_correlation = true;
00095 }
00096 else if (sarg=="-make_plots") {
00097 cout << "Will make files suitable for plotting\n"<<endl;
00098 make_plots = true;
00099 }
00100 else if (sarg.substr(0,1)=="-") modeList.push_back(sarg);
00101 else paramFile = sarg;
00102 }
00103
00104
00105 if (paramFile=="") {
00106 MESSAGE<<"No inputFile parameter file given!!"<<ENDM_FATAL;
00107 getHelp();
00108 }
00109 parameterMap inMap;
00110 ReadParsFromFile(inMap, paramFile);
00111
00112
00113 string oscarFile = inMap.getItem( "oscar_file", (string)"thermal_gauss01.dat" );
00114 vector<COSCARLine> particleList = filterPID(
00115 readOSCARFile(oscarFile),
00116 inMap.getItem( "pid", 211 )
00117 );
00118
00119
00120
00121 for (vector<string>::iterator mode=modeList.begin(); mode!=modeList.end(); ++mode) {
00122 if (*mode == "-1d") make1dSourceAndCorrelation( particleList, inMap, skip_correlation, make_plots );
00123 else if (*mode == "-1doldway") make1dSourceAndCorrelationOldWay( particleList, inMap, skip_correlation, make_plots );
00124 else if (*mode == "-3dcart") make3dCartSourceAndCorrelation( particleList, inMap, skip_correlation, make_plots );
00125 else if (*mode == "-3dcart_ex") make3dCartExpandedSourceAndCorrelation( particleList, inMap, skip_correlation, make_plots );
00126 else if (*mode == "-3dsphr") make3dSphrSourceAndCorrelation( particleList, inMap, skip_correlation, make_plots );
00127 else if (*mode == "-crab" && !skip_correlation) makeCRABCorrelations( particleList, inMap, make_plots );
00128 else if (*mode == "-sourcemaker") makeSourceMakerSources( particleList, inMap, make_plots );
00129 else MESSAGE<< "Unknown mode: "<<*mode<<ENDM_WARN;
00130 cout << endl;
00131 }
00132 return true;
00133 }
00134
00135
00136
00137
00138
00139
00140 void make1dSourceAndCorrelation( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots ){
00141 cout << "\nGenerating 1d source"<<endl;
00142 CSourceFtn1dHisto souin = getSource1d( particleList, inMap );
00143 cout << endl;
00144 parameterMap coutMap, soutMap;
00145 souin.Write(soutMap);
00146 WriteParsToFile(soutMap,parameter::getS(inMap, "source_file", "output_source_1d.dat"));
00147 if (make_plots){
00148 plot1d( souin, getS( inMap, "source_plot_filename", "output_source_1d_plot.dat" ));
00149 }
00150 if (!skip_correlation) {
00151 cout << "\nGenerating 1d correlation from source"<<endl;
00152 CCorrFtn1dHisto answer = getCorrelation1d( particleList, inMap );
00153 cout << endl;
00154 answer.Write(coutMap);
00155 WriteParsToFile(coutMap,parameter::getS(inMap, "correlation_file", "output_correlation_1d.dat"));
00156 if (make_plots){
00157 plot1d( answer, getS( inMap, "correlation_plot_filename", "output_correlation_1d_plot.dat" ));
00158 }
00159 }
00160 }
00161
00162
00163 void make1dSourceAndCorrelationOldWay( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots ){
00164 cout << "\nGenerating 1d source, the old way"<<endl;
00165 parameterMap souMapIn;
00166 souMapIn = parameter::getMap(inMap, "source_settings");
00167 CSourceFtn1dHisto soualt = getSource1dOldWay( particleList, souMapIn );
00168 cout << endl;
00169 parameterMap coutMap, soutMap;
00170 soualt.Write(soutMap);
00171 WriteParsToFile(soutMap,parameter::getS(inMap, "source_file", "output_source_1d.dat"));
00172 if (make_plots){
00173 plot1d( soualt, getS( inMap, "source_plot_filename", "output_source_1d_plot.dat" ));
00174 }
00175 if (!skip_correlation) {
00176 cout << "\nGenerating 1d correlation from source"<<endl;
00177 CCorrFtn1dHisto answer;
00178 CBasisFuncImager1d imager;
00179 imager.convertSourceToCorrelation(soualt, answer, inMap);
00180 cout << endl;
00181 answer.Write(coutMap);
00182 WriteParsToFile(coutMap,parameter::getS(inMap, "correlation_file", "output_correlation_1d.dat"));
00183 if (make_plots){
00184 plot1d( answer, getS( inMap, "correlation_plot_filename", "output_correlation_1d_plot.dat" ));
00185 }
00186 }
00187 }
00188
00189
00190 void make3dCartSourceAndCorrelation( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots ){
00191 cout << "\nGenerating 3d source, in Cartesian"<<endl;
00192 CSourceFtn3dHisto sanswer = getSource3dCart( particleList, inMap );
00193 cout << endl;
00194 parameterMap coutMap, soutMap;
00195 sanswer.Write(soutMap);
00196 WriteParsToFile( soutMap, parameter::getS(inMap, "source_file", "output_source_3d.dat") );
00197 if (make_plots){
00198 plot3dSlices( sanswer,
00199 getS( inMap, "source_3d_side_plot_filename", "output_source_3d_side_plot.dat" ),
00200 getS( inMap, "source_3d_out_plot_filename", "output_source_3d_out_plot.dat" ),
00201 getS( inMap, "source_3d_long_plot_filename", "output_source_3d_long_plot.dat" )
00202 );
00203 }
00204 if (!skip_correlation) {
00205 cout << "\nGenerating 3d correlation, in Cartesian"<<endl;
00206 CCorrFtn3dHisto answer = getCorrelation3dCart( particleList, inMap );
00207 cout << endl;
00208 answer.Write(coutMap);
00209 WriteParsToFile( coutMap, parameter::getS(inMap, "correlation_file", "output_correlation_3d.dat") );
00210 if (make_plots){
00211 plot3dSlices( answer,
00212 getS( inMap, "correlation_3d_side_plot_filename", "output_correlation_3d_side_plot.dat" ),
00213 getS( inMap, "correlation_3d_out_plot_filename", "output_correlation_3d_out_plot.dat" ),
00214 getS( inMap, "correlation_3d_long_plot_filename", "output_correlation_3d_long_plot.dat" )
00215 );
00216 }
00217 }
00218 }
00219
00220
00221 void make3dCartExpandedSourceAndCorrelation( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots ){
00222 cout << "\nGenerating 3d source in Cartesian coordinates"<<endl;
00223 CSourceFtn3dHisto souin3d = getSource3dCart( particleList, inMap );
00224 cout << "\nExpanding 3d source into Ylm's"<<endl;
00225 CSourceFtn3dSphr< CSourceFtn1dHisto > soutemp = expand( souin3d, inMap );
00226 cout << endl;
00227 parameterMap coutMap, soutMap, souinMap;
00228 souin3d.Write(souinMap);
00229 soutemp.Write(soutMap);
00230 WriteParsToFile( souinMap, parameter::getS(inMap, "source_file_cart", "output_source_3d_cart.dat") );
00231 WriteParsToFile( soutMap, parameter::getS(inMap, "source_file", "output_source_3d_sphr.dat") );
00232 soutemp.writeTerms();
00233 if (make_plots){
00234 cout << "Do nothing for plotting"<<endl;
00235 }
00236 if (!skip_correlation) {
00237 cout << "\nGenerating 3d correlation from source"<<endl;
00238 CCorrFtn3dSphr answer;
00239 UncoupledHistoImager3d imager;
00240 imager.convertSourceToCorrelation(soutemp, answer, inMap);
00241 cout << endl;
00242 answer.Write(coutMap);
00243 WriteParsToFile( coutMap, parameter::getS(inMap, "correlation_file", "output_correlation_3d_cart.dat") );
00244 answer.writeTerms( );
00245 if (make_plots){
00246 cout << "Do nothing for plotting"<<endl;
00247 }
00248 }
00249 }
00250
00251
00252 void make3dSphrSourceAndCorrelation( vector<COSCARLine> particleList, parameterMap& inMap, bool skip_correlation, bool make_plots ){
00253 cout << "\nGenerating 3d source directly in Ylm's"<<endl;
00254 CSourceFtn3dSphr< CSourceFtn1dHisto > soualt3d = getSource3dSphr( particleList, inMap );
00255 cout << endl;
00256 parameterMap soutMap, coutMap;
00257 soualt3d.Write(soutMap);
00258 soualt3d.writeTerms();
00259 WriteParsToFile( soutMap, getS( inMap, "source_file", "output_source_3dsphr.dat" ) );
00260 if (make_plots){
00261 plot3dSlices( soualt3d,
00262 getS( inMap, "source_3d_side_plot_filename", "output_source_3d_side_plot.dat" ),
00263 getS( inMap, "source_3d_out_plot_filename", "output_source_3d_out_plot.dat" ),
00264 getS( inMap, "source_3d_long_plot_filename", "output_source_3d_long_plot.dat" )
00265 );
00266 }
00267 if (!skip_correlation) {
00268 cout << "\nGenerating 3d correlation from source"<<endl;
00269 CCorrFtn3dSphr answer;
00270 UncoupledHistoImager3d imager;
00271 imager.convertSourceToCorrelation(soualt3d, answer, inMap);
00272 cout << endl;
00273 answer.Write(coutMap);
00274 answer.writeTerms();
00275 WriteParsToFile( soutMap, getS( inMap, "correlation_file", "output_correlation_3dsphr.dat" ) );
00276 if (make_plots){
00277 plot3dSlices( answer,
00278 getS( inMap, "correlation_3d_side_plot_filename", "output_correlation_3d_side_plot.dat" ),
00279 getS( inMap, "correlation_3d_out_plot_filename", "output_correlation_3d_out_plot.dat" ),
00280 getS( inMap, "correlation_3d_long_plot_filename", "output_correlation_3d_long_plot.dat" )
00281 );
00282 }
00283 }
00284 }
00285
00286
00287 void makeCRABCorrelations( vector<COSCARLine> particleList, parameterMap& inMap, bool make_plots ){
00288 cout << "\nGenerating 1d correlation like CRAB used to do"<<endl;
00289 COSCARCorrelationGenerator1d corr1dgen;
00290 CCorrFtn1dHisto corr1d(corr1dgen.generateCorrelation(particleList, inMap));
00291 cout << endl;
00292 cout << "\nGenerating 3d correlation like CRAB used to do"<<endl;
00293 COSCARCorrelationGenerator3dCart corr3dgen;
00294 CCorrFtn3dHisto corr3d(corr3dgen.generateCorrelation(particleList, inMap));
00295 cout << endl;
00296 TNT::Array1D< int > numPairs1d = corr1dgen.pairCount;
00297 TNT::Array1D< int > numPairs3d = corr3dgen.pairCount;
00298 cout << "\nOutputting 1d and 3d correlations like CRAB used to do"<<endl;
00299 ofstream outFile1d(parameter::getS(inMap, "crab_correlation_1d_filename", "correlation_qinv.dat").c_str());
00300 outFile1d << "!As a function of reduced momentum, k=(p2-p1)/2:\n"
00301 << "!k(MeV/c) C(k) +/- numcounts&dencounts"<<endl;
00302 for (int i=0; i<corr1d.ndata; ++i)
00303 {
00304 outFile1d<< setw(6) << setprecision(3) << corr1d.midBin(i) << " "
00305 << setw(6) << setprecision(3) << corr1d.data[i] << " "
00306 << setw(9) << setprecision(6) << corr1d.uncert[i] << " "
00307 << numPairs1d[i] << " " << 0 << endl;
00308 }
00309 ofstream outFile3d(parameter::getS(inMap, "crab_correlation_3d_filename", "correlation_qinv3d.dat").c_str());
00310 outFile3d << "!As a function of reduced momentum, k=(p2-p1)/2:\n"
00311 << "! kx ky kz C(k) +/- numcounts&dencounts"<<endl;
00312 for (int ix=0; ix<corr3d.nx; ++ix)
00313 for (int iy=0; iy<corr3d.ny; ++iy)
00314 for (int iz=0; iz<corr3d.nz; ++iz)
00315 {
00316 int i = corr3d.whatIndex(ix,iy,iz);
00317 outFile3d
00318 << setw(6) << setprecision(3) << corr3d.midBinX(ix) << " "
00319 << setw(6) << setprecision(3) << corr3d.midBinY(iy) << " "
00320 << setw(6) << setprecision(3) << corr3d.midBinZ(iz) << " "
00321 << setw(6) << setprecision(3) << corr3d.data[i] << " "
00322 << setw(9) << setprecision(6) << corr3d.uncert[i] << " "
00323 << numPairs3d[i] << " " << 0 << endl;
00324 }
00325 if (make_plots) {
00326 plot1d( corr1d, getS( inMap, "correlation_plot_filename", "output_correlation_1d_plot.dat" ) );
00327 plot3dSlices( corr3d,
00328 getS( inMap, "correlation_3d_side_plot_filename", "output_correlation_3d_side_plot.dat" ),
00329 getS( inMap, "correlation_3d_out_plot_filename", "output_correlation_3d_out_plot.dat" ),
00330 getS( inMap, "correlation_3d_long_plot_filename", "output_correlation_3d_long_plot.dat" )
00331 );
00332 }
00333 }
00334
00335
00336 void makeSourceMakerSources( vector<COSCARLine> particleList, parameterMap& inMap, bool make_plots ){
00337 cout << "\nGenerating 1d source slices like SourceMaker used to do"<<endl;
00338 CSourceFtn3dHisto source3d = getSource3dCart( particleList, inMap );
00339 cout << endl;
00340 cout << "\nGenerating 3d source slices like SourceMaker used to do"<<endl;
00341 CSourceFtn1dHisto source1d = getSource1d( particleList, inMap );
00342 if (make_plots) {
00343 plot1d( source1d, getS( inMap, "source_1d_plot_filename", "output_source_1d_plot.dat" ));
00344 cout << endl;
00345 plot3dSlices( source3d,
00346 getS( inMap, "source_3d_side_plot_filename", "output_source_3d_side_plot.dat" ),
00347 getS( inMap, "source_3d_out_plot_filename", "output_source_3d_out_plot.dat" ),
00348 getS( inMap, "source_3d_long_plot_filename", "output_source_3d_long_plot.dat" ) );
00349 }
00350 }
00351
00352
00354 CSourceFtn1dHisto getSource1dOldWay( vector<COSCARLine> particleList, parameterMap m ){
00355 CSourceFtn1dHisto souin(getOSCARSource1d(particleList, m));
00356 return souin;
00357 }
00358
00359
00361 CSourceFtn1dHisto getSource1d( vector<COSCARLine> particleList, parameterMap m ){
00362 COSCARSourceGenerator1d sougen;
00363 CSourceFtn1dHisto soualt(sougen.generateSource(particleList, m));
00364 return soualt;
00365 }
00366
00367
00369 CCorrFtn1dHisto getCorrelation1d( vector<COSCARLine> particleList, parameterMap m ){
00370 COSCARCorrelationGenerator1d corrgen;
00371 CCorrFtn1dHisto corralt(corrgen.generateCorrelation(particleList, m));
00372 return corralt;
00373 }
00374
00375
00377 CCorrFtn1dHisto getCorrelation1dOldWay( const CSourceFtn1dHisto& souin, parameterMap m ){
00378 CBasisFuncImager1d imager;
00379 CCorrFtn1dHisto answer;
00380 imager.convertSourceToCorrelation( souin, answer, m );
00381 return answer;
00382 }
00383
00384
00385 CSourceFtn3dSphr< CSourceFtn1dHisto > getSource3dSphr( vector<COSCARLine> particleList, parameterMap m ){
00386 COSCARSourceGenerator3dSphr sougen;
00387 CSourceFtn3dSphr< CSourceFtn1dHisto > souin(sougen.generateSource(particleList, m));
00388 return souin;
00389 }
00390
00391
00392 CSourceFtn3dHisto getSource3dCart( vector<COSCARLine> particleList, parameterMap m ){
00393 COSCARSourceGenerator3dCart sougen;
00394 CSourceFtn3dHisto souin(sougen.generateSource(particleList, m));
00395 return souin;
00396 }
00397
00398
00399 CCorrFtn3dHisto getCorrelation3dCart( vector<COSCARLine> particleList, parameterMap m ){
00400 COSCARCorrelationGenerator3dCart sougen;
00401 CCorrFtn3dHisto souin(sougen.generateCorrelation(particleList, m));
00402 return souin;
00403 }
00404
00405
00406 void plot1d( CObject1d& obj, string outFileName )
00407 {
00408 cout << "Writing 1d plot file "<<outFileName<<endl;
00409 ofstream dataOut( outFileName.c_str() );
00410 for (int i=0; i<52; ++i){
00411 double x = (double)i;
00412 dataOut << x << " " << obj.getValue(x) << " " << obj.getError(x)<<endl;
00413 }
00414 }
00415
00416
00417 void plot1d( CHistogram1d& obj, string outFileName )
00418 {
00419 cout << "Writing 1d plot file "<<outFileName<<endl;
00420 ofstream dataOut( outFileName.c_str() );
00421 for (int i=0; i<obj.ndata; ++i){
00422 double x = obj.midBin(i);
00423 dataOut << x << " " << obj.getValue(x) << " " << obj.getError(x)<<endl;
00424 }
00425 }
00426
00427
00428 void plot3dSlices( CObject3d& obj, string xFileName, string yFileName, string zFileName ){
00429 cout << "Writing X slice to plot file "<<xFileName<<endl;
00430 ofstream dataX( xFileName.c_str() );
00431 for (int i=0; i<52; ++i){
00432 double x = (double)i;
00433 dataX << x << " " << obj.getValueCart(x,0.,0.) << " " << obj.getErrorCart(x,0.,0.)<<endl;
00434 }
00435 cout << "Writing Y slice to plot file "<<yFileName<<endl;
00436 ofstream dataY( yFileName.c_str() );
00437 for (int i=0; i<52; ++i){
00438 double x = (double)i;
00439 dataY << x << " " << obj.getValueCart(0.,x,0.) << " " << obj.getErrorCart(0.,x,0.)<<endl;
00440 }
00441 cout << "Writing Z slice to plot file "<<zFileName<<endl;
00442 ofstream dataZ( zFileName.c_str() );
00443 for (int i=0; i<52; ++i){
00444 double x = (double)i;
00445 dataZ << x << " " << obj.getValueCart(0.,0.,x) << " " << obj.getErrorCart(0.,0.,x)<<endl;
00446 }
00447 }
00448
00449
00450 void plot3dSlices( CHistogram3d& obj, string xFileName, string yFileName, string zFileName ){
00451 cout << "Writing X slice to plot file "<<xFileName<<endl;
00452 ofstream dataX( xFileName.c_str() );
00453 for (int i=0; i<obj.nx; ++i){
00454 double x = obj.midBinX(i);
00455 dataX << x << " " << obj.getValueCart(x,0.,0.) << " " << obj.getErrorCart(x,0.,0.)<<endl;
00456 }
00457 cout << "Writing Y slice to plot file "<<yFileName<<endl;
00458 ofstream dataY( yFileName.c_str() );
00459 for (int i=0; i<obj.ny; ++i){
00460 double x = obj.midBinY(i);
00461 dataY << x << " " << obj.getValueCart(0.,x,0.) << " " << obj.getErrorCart(0.,x,0.)<<endl;
00462 }
00463 cout << "Writing Z slice to plot file "<<zFileName<<endl;
00464 ofstream dataZ( zFileName.c_str() );
00465 for (int i=0; i<obj.nz; ++i){
00466 double x = obj.midBinZ(i);
00467 dataZ << x << " " << obj.getValueCart(0.,0.,x) << " " << obj.getErrorCart(0.,0.,x)<<endl;
00468 }
00469 }