diff --git a/AnalysisLibrary.h b/AnalysisLibrary.h new file mode 100644 index 0000000..baf2118 --- /dev/null +++ b/AnalysisLibrary.h @@ -0,0 +1,235 @@ +#ifndef Analysis_Library_h +#define Analysis_Library_h + +#include +#include +#include +#include +#include +#include +#include + +std::vector SplitStr(std::string tempLine, std::string splitter, int shift = 0){ + + std::vector output; + + size_t pos; + do{ + pos = tempLine.find(splitter); // fine splitter + if( pos == 0 ){ //check if it is splitter again + tempLine = tempLine.substr(pos+1); + continue; + } + + std::string secStr; + if( pos == std::string::npos ){ + secStr = tempLine; + }else{ + secStr = tempLine.substr(0, pos+shift); + tempLine = tempLine.substr(pos+shift); + } + + //check if secStr is begin with space + while( secStr.substr(0, 1) == " "){ + secStr = secStr.substr(1); + }; + + //check if secStr is end with space + while( secStr.back() == ' '){ + secStr = secStr.substr(0, secStr.size()-1); + } + + output.push_back(secStr); + //printf(" |%s---\n", secStr.c_str()); + + }while(pos != std::string::npos ); + + return output; +} + +std::vector> combination(std::vector arr, int r){ + + std::vector> output; + + int n = arr.size(); + std::vector v(n); + std::fill(v.begin(), v.begin()+r, 1); + do { + //for( int i = 0; i < n; i++) { printf("%d ", v[i]); }; printf("\n"); + + std::vector temp; + for (int i = 0; i < n; ++i) { + if (v[i]) { + //printf("%.1f, ", arr[i]); + temp.push_back(arr[i]); + } + } + //printf("\n"); + + output.push_back(temp); + + } while (std::prev_permutation(v.begin(), v.end())); + + return output; +} + +double* sumMeanVar(std::vector data){ + + int n = data.size(); + double sum = 0; + for( int k = 0; k < n; k++) sum += data[k]; + double mean = sum/n; + double var = 0; + for( int k = 0; k < n; k++) var += pow(data[k] - mean,2); + + static double output[3]; + output[0] = sum; + output[1] = mean; + output[2] = var; + + return output; +} + +double* fitSlopeIntercept(std::vector dataX, std::vector dataY){ + + double * smvY = sumMeanVar(dataY); + double sumY = smvY[0]; + double meanY = smvY[1]; + + double * smvX = sumMeanVar(dataX); + double sumX = smvX[0]; + double meanX = smvX[1]; + double varX = smvX[2]; + + int n = dataX.size(); + double sumXY = 0; + for( int j = 0; j < n; j++) sumXY += dataX[j] * dataY[j]; + + double slope = ( sumXY - sumX * sumY/n ) / varX; + double intercept = meanY - slope * meanX; + + static double output[2]; + output[0] = slope; + output[1] = intercept; + + return output; + +} + + + +std::vector> FindMatchingPair(std::vector enX, std::vector enY){ + + //output[0] = fitEnergy; + //output[1] = refEnergy; + + int nX = enX.size(); + int nY = enY.size(); + + std::vector fitEnergy; + std::vector refEnergy; + + if( nX > nY ){ + + std::vector> output = combination(enX, nY); + + double * smvY = sumMeanVar(enY); + double sumY = smvY[0]; + double meanY = smvY[1]; + double varY = smvY[2]; + + double optRSquared = 0; + double absRSqMinusOne = 1; + int maxID = 0; + + for( int k = 0; k < (int) output.size(); k++){ + + double * smvX = sumMeanVar(output[k]); + double sumX = smvX[0]; + double meanX = smvX[1]; + double varX = smvX[2]; + + double sumXY = 0; + for( int j = 0; j < nY; j++) sumXY += output[k][j] * enY[j]; + + double rSq = abs(sumXY - sumX*sumY/nY)/sqrt(varX*varY); + + //for( int j = 0; j < nY ; j++){ printf("%.1f, ", output[k][j]); }; printf("| %.10f\n", rSq); + + if( abs(rSq-1) < absRSqMinusOne ) { + absRSqMinusOne = abs(rSq-1); + optRSquared = rSq; + maxID = k; + } + } + + fitEnergy = output[maxID]; + refEnergy = enY; + + printf(" R^2 : %.20f\n", optRSquared); + + //calculation fitting coefficient + //double * si = fitSlopeIntercept(fitEnergy, refEnergy); + //printf( " y = %.4f x + %.4f\n", si[0], si[1]); + + }else if( nX < nY ){ + + std::vector> output = combination(enY, nX); + + + double * smvX = sumMeanVar(enX); + double sumX = smvX[0]; + double meanX = smvX[1]; + double varX = smvX[2]; + + double optRSquared = 0; + double absRSqMinusOne = 1; + int maxID = 0; + + for( int k = 0; k < (int) output.size(); k++){ + + double * smvY = sumMeanVar(output[k]); + double sumY = smvY[0]; + double meanY = smvY[1]; + double varY = smvY[2]; + + double sumXY = 0; + for( int j = 0; j < nX; j++) sumXY += output[k][j] * enX[j]; + + double rSq = abs(sumXY - sumX*sumY/nX)/sqrt(varX*varY); + + //for( int j = 0; j < nX ; j++){ printf("%.1f, ", output[k][j]); }; printf("| %.10f\n", rSq); + + if( abs(rSq-1) < absRSqMinusOne ) { + absRSqMinusOne = abs(rSq-1); + optRSquared = rSq; + maxID = k; + } + } + + fitEnergy = enX; + refEnergy = output[maxID]; + printf(" R^2 : %.20f\n", optRSquared); + + }else{ + fitEnergy = enX; + refEnergy = enY; + + //if nX == nY, ther could be cases that only partial enX and enY are matched. + + } + + + printf("fitEnergy = ");for( int k = 0; k < (int) fitEnergy.size() ; k++){ printf("%7.2f, ", fitEnergy[k]); }; printf("\n"); + printf("refEnergy = ");for( int k = 0; k < (int) refEnergy.size() ; k++){ printf("%7.2f, ", refEnergy[k]); }; printf("\n"); + + std::vector> haha; + haha.push_back(fitEnergy); + haha.push_back(refEnergy); + + return haha; + +} + +#endif + diff --git a/Analyzer.C b/Analyzer.C index 785dfec..5b776ed 100644 --- a/Analyzer.C +++ b/Analyzer.C @@ -11,12 +11,12 @@ //############################################ User setting -int rawEnergyRange[2] = {100, 6000}; // in ch +int rawEnergyRange[2] = {500, 6000}; // in ch int energyRange[3] = {1, 100, 2000}; // keV {resol, min, max} double BGO_threshold = 100; // in ch -TString e_corr = "e_corr.txt"; +TString e_corr = "correction_e.dat"; //############################################ end of user setting @@ -30,15 +30,15 @@ vector> eCorr; //############################################ histogram declaration TH2F * heVID; -TH1F * he[NCLOVER]; +TH1F * he[NCRYSTAL]; -TH2F * hgg[NCLOVER][NCLOVER]; +TH2F * hgg[NCRYSTAL][NCRYSTAL]; TH2F * hcoin; ///----- after calibration and BGO veto TH2F * heCalVID; -TH1F * heCal[NCLOVER]; +TH1F * heCal[NCRYSTAL]; TH2F * hcoinBGO; void Analyzer::Begin(TTree * tree){ @@ -49,23 +49,23 @@ void Analyzer::Begin(TTree * tree){ printf("======================== Histograms declaration\n"); - heVID = new TH2F("heVID", "e vs ID; det ID; e [ch]", NCLOVER, 0, NCLOVER, rawEnergyRange[1] - rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); - heCalVID = new TH2F("heCalVID", Form("eCal vs ID (BGO veto > %.1f); det ID; e [ch]", BGO_threshold), NCLOVER, 0, NCLOVER, rawEnergyRange[1] - rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); - for( int i = 0; i < NCLOVER; i ++){ + heVID = new TH2F("heVID", "e vs ID; det ID; e [ch]", NCRYSTAL, 0, NCRYSTAL, rawEnergyRange[1] - rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); + heCalVID = new TH2F("heCalVID", Form("eCal vs ID (BGO veto > %.1f); det ID; Energy [keV]", BGO_threshold), NCRYSTAL, 0, NCRYSTAL, (energyRange[2] - energyRange[1])/energyRange[0], energyRange[1], energyRange[2]); + for( int i = 0; i < NCRYSTAL; i ++){ he[i] = new TH1F( Form("he%02d", i), Form("e -%02d", i), rawEnergyRange[1] - rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); - heCal[i] = new TH1F(Form("heCal%02d", i), Form("e -%02d (BGO veto > %.1f)", i, BGO_threshold), rawEnergyRange[1] - rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); + heCal[i] = new TH1F(Form("heCal%02d", i), Form("eCal -%02d (BGO veto > %.1f); Energy [keV]; count / %d keV", i, BGO_threshold, energyRange[0]), (energyRange[2] - energyRange[1])/energyRange[0], energyRange[1], energyRange[2]); } - for( int i = 0; i < NCLOVER; i++){ - for( int j = i+1; j < NCLOVER; j++){ + for( int i = 0; i < NCRYSTAL; i++){ + for( int j = i+1; j < NCRYSTAL; j++){ //hgg[i][j] = new TH2F(Form("hgg%02d%02d", i, j), Form("e%02d vs e%02d; e%02d; e%02d", i, j, i, j), // (rawEnergyRange[1] - rawEnergyRange[0])/2, rawEnergyRange[0], rawEnergyRange[1], // (rawEnergyRange[1] - rawEnergyRange[0])/2, rawEnergyRange[0], rawEnergyRange[1]); } } - hcoin = new TH2F("hcoin", "detector coin.; det ID; det ID", NCLOVER, 0, NCLOVER, NCLOVER, 0 , NCLOVER); - hcoinBGO = new TH2F("hcoinBGO", Form("detector coin. (BGO veto > %.1f); det ID; det ID", BGO_threshold), NCLOVER, 0, NCLOVER, NCLOVER, 0 , NCLOVER); + hcoin = new TH2F("hcoin", "detector coin.; det ID; det ID", NCRYSTAL, 0, NCRYSTAL, NCRYSTAL, 0 , NCRYSTAL); + hcoinBGO = new TH2F("hcoinBGO", Form("detector coin. (BGO veto > %.1f); det ID; det ID", BGO_threshold), NCRYSTAL, 0, NCRYSTAL, NCRYSTAL, 0 , NCRYSTAL); printf("======================== End of histograms declaration\n"); @@ -102,7 +102,7 @@ Bool_t Analyzer::Process(Long64_t entry){ if( multi == 0 ) return kTRUE; ///=========== Looping Crystals - for( int detID = 0; detID < NCLOVER ; detID ++){ + for( int detID = 0; detID < NCRYSTAL ; detID ++){ //======== baics gate when no energy or pileup if( TMath::IsNaN(e[detID])) continue; @@ -113,7 +113,7 @@ Bool_t Analyzer::Process(Long64_t entry){ he[detID]->Fill(e[detID]); - for( int detJ = detID +1; detJ < NCLOVER; detJ++) { + for( int detJ = detID +1; detJ < NCRYSTAL; detJ++) { if( TMath::IsNaN(e[detJ])) continue; //hgg[detID][detJ]->Fill(e[detID], e[detJ]); // x then y hcoin->Fill(detID, detJ); @@ -139,7 +139,7 @@ Bool_t Analyzer::Process(Long64_t entry){ heCalVID->Fill( detID, eCal); heCal[detID]->Fill(eCal); - for( int detJ = detID +1; detJ < NCLOVER; detJ++) { + for( int detJ = detID +1; detJ < NCRYSTAL; detJ++) { if( TMath::IsNaN(e[detJ])) continue; hcoinBGO->Fill(detID, detJ); } @@ -177,13 +177,14 @@ void Analyzer::Terminate(){ cCanvas->cd(4); cCanvas->cd(4)->SetLogz(1); hcoinBGO->Draw("colz"); - + + printf("=============== loaded AutoFit.C, try showFitMethos()\n"); + gROOT->ProcessLine(".L AutoFit.C"); printf("=============== Analyzer Utility\n"); gROOT->ProcessLine(".L Analyzer_Utili.c"); gROOT->ProcessLine("listDraws()"); - printf("=============== loaded AutoFit.C\n"); - gROOT->ProcessLine(".L AutoFit.C"); + } diff --git a/Analyzer.h b/Analyzer.h index 96f24c1..6d3d8e3 100644 --- a/Analyzer.h +++ b/Analyzer.h @@ -15,6 +15,7 @@ #include #include "mapping.h" +#include "AnalysisLibrary.h" // Header file for the classes stored in the TTree if any. @@ -26,9 +27,9 @@ public : // Declaration of leaf types ULong64_t evID; - Double_t e[NCLOVER]; - ULong64_t t[NCLOVER]; - UShort_t p[NCLOVER]; + Double_t e[NCRYSTAL]; + ULong64_t t[NCRYSTAL]; + UShort_t p[NCRYSTAL]; Double_t bgo[NBGO]; Double_t other[NOTHER]; Int_t multi; @@ -112,49 +113,9 @@ void Analyzer::SlaveTerminate(){ } - - -std::vector SplitStr(std::string tempLine, std::string splitter, int shift = 0){ - - std::vector output; - - size_t pos; - do{ - pos = tempLine.find(splitter); /// fine splitter - if( pos == 0 ){ ///check if it is splitter again - tempLine = tempLine.substr(pos+1); - continue; - } - - std::string secStr; - if( pos == std::string::npos ){ - secStr = tempLine; - }else{ - secStr = tempLine.substr(0, pos+shift); - tempLine = tempLine.substr(pos+shift); - } - - ///check if secStr is begin with space - while( secStr.substr(0, 1) == " "){ - secStr = secStr.substr(1); - }; - - ///check if secStr is end with space - while( secStr.back() == ' '){ - secStr = secStr.substr(0, secStr.size()-1); - } - - output.push_back(secStr); - //printf(" |%s---\n", secStr.c_str()); - - }while(pos != std::string::npos ); - - return output; -} - std::vector> LoadCorrectionParameters(TString corrFile){ - printf("==================== load correction parameters : %s", corrFile.Data()); + printf(" load correction parameters : %s", corrFile.Data()); std::ifstream file; file.open(corrFile.Data()); @@ -192,13 +153,17 @@ std::vector> LoadCorrectionParameters(TString corrFile){ printf("det : %2d | ", i ); int len = (int) corr[i].size(); for( int j = 0; j < len - 1 ; j++){ - printf("%6.2f, ", corr[i][j]); + printf("%14.6f, ", corr[i][j]); } - printf("%6.2f\n", corr[i][len-1]); + printf("%14.6f\n", corr[i][len-1]); } }else{ - printf(".... fail\n"); + printf(".... fail\n"); + std::vector temp = {0, 1}; + for( int i = 0; i < NCRYSTAL; i++){ + corr.push_back(temp); + } } return corr; diff --git a/Analyzer_Utili.c b/Analyzer_Utili.c index 7823316..5a8f93d 100644 --- a/Analyzer_Utili.c +++ b/Analyzer_Utili.c @@ -2,9 +2,11 @@ void listDraws(void) { printf("------------------- List of Plots -------------------\n"); printf(" newCanvas() - Create a new Canvas\n"); printf("-----------------------------------------------------\n"); - printf(" rawEvID() - Raw e vs ID\n"); - printf(" drawE() - Raw e for all %d detectors\n", NCLOVER); - //printf(" drawGG() - Gamma - Gamma Coincident for all %d detectors\n", NCLOVER); + printf(" rawEvID() - e vs ID\n"); + printf(" drawE() - e for all %d detectors\n", NCRYSTAL); + //printf(" drawGG() - Gamma - Gamma Coincident for all %d detectors\n", NCRYSTAL); + printf("-----------------------------------------------------\n"); + printf(" energyCalibration() - Calibrate energy \n"); printf("-----------------------------------------------------\n"); } @@ -16,58 +18,85 @@ void newCanvas(int sizeX = 800, int sizeY = 600, int posX = 0, int posY = 0){ cNewCanvas->cd(); } -void rawEvID(){ +void rawEvID(bool cal = false){ TCanvas * cRawID = (TCanvas *) gROOT->FindObjectAny("cRawID"); if( cRawID == NULL ) cRawID = new TCanvas("cRawID", "raw ID", 1000, 800); cRawID->cd(1)->SetGrid(); - heVID->Draw("colz"); + cal ? heCalVID->Draw("colz") : heVID->Draw("colz"); } -void drawE(bool isLogy = false, bool cali = false){ +void drawE(int CloverID = -1, bool cali = false, bool isLogy = false, double xMin = 0, double xMax = 0){ - int nCrystal = 4; - int numCol = NCLOVER / nCrystal; + int nCrystalPerClover = 4; + int nClover = NCRYSTAL / nCrystalPerClover; + + if( CloverID >= nClover ) { + printf("Clover-ID > nClover = %d. \n", nClover); + return; + } int size = 300; TCanvas *cRawE = (TCanvas *) gROOT->FindObjectAny("cRawE"); - if( cRawE == NULL ) cRawE = new TCanvas("cRawE", cali ? "Cal e" : "Raw e", size * numCol, size * nCrystal); - cRawE->Clear();cRawE->Divide(numCol, 4, 0); + if( cRawE == NULL ) cRawE = new TCanvas("cRawE", cali ? "Cal e" : "Raw e", size * nClover, size * nCrystalPerClover); + cRawE->Clear(); + if( CloverID >= 0 ) nClover = 1; + cRawE->Divide(nClover, nCrystalPerClover, 0); - for (Int_t i = 0; i < nCrystal; i++) { - for( Int_t j = 0; j < numCol; j++){ - int canvasID = numCol * i + j + 1; + + ///find max y + double maxY = 0; + int nDet = nClover*nCrystalPerClover; + for( int i = (CloverID < 0 ? 0 : nCrystalPerClover*CloverID) ; i < (CloverID < 0 ? nDet : nCrystalPerClover*CloverID + nDet) ; i++){ + int mBin = cali ? heCal[i]->GetMaximumBin() : he[i]->GetMaximumBin(); + double max = cali ? heCal[i]->GetBinContent(mBin) : he[i]->GetBinContent(mBin); + if( max > maxY ) maxY = max; + } + maxY = maxY * 1.1; + ///printf("max Y : %f \n", maxY); + + for (Int_t i = 0; i < nClover; i++) { + for( Int_t j = 0; j < nCrystalPerClover; j++){ + int canvasID = CloverID < 0 ? nClover*j+ i + 1 : j + 1; cRawE->cd(canvasID); cRawE->cd(canvasID)->SetGrid(); cRawE->cd(canvasID)->SetTickx(2); cRawE->cd(canvasID)->SetTicky(2); + cRawE->cd(canvasID)->SetBottomMargin(0.06); + if ( i == nClover -1 ) cRawE->cd(canvasID)->SetRightMargin(0.002); if( isLogy ) cRawE->cd(canvasID)->SetLogy(); - int hID = nCrystal*j+ i; + int hID = CloverID < 0 ? nCrystalPerClover*i+ j : nCrystalPerClover * CloverID + j ; if( cali ) { + if ( xMin != 0 || xMax != 0 ) heCal[hID]->GetXaxis()->SetRangeUser(xMin, xMax); + heCal[hID]->SetMaximum(maxY); heCal[hID]->Draw(""); }else{ + if ( xMin != 0 || xMax != 0 ) he[hID]->GetXaxis()->SetRangeUser(xMin, xMax); + he[hID]->SetMaximum(maxY); he[hID]->Draw(""); } } } - + + cRawE->SetCrosshair(1); + } /** void drawGG(){ int nCrystal = 4; - int numCol = NCLOVER / nCrystal; + int numCol = NCRYSTAL / nCrystal; int size = 300; TCanvas *cGG = (TCanvas *) gROOT->FindObjectAny("cGG"); - if( cGG == NULL ) cGG = new TCanvas("cGG", "Gamma - Gamma Coin.", size * NCLOVER, size * NCLOVER); - cGG->Clear();cGG->Divide(NCLOVER, NCLOVER); + if( cGG == NULL ) cGG = new TCanvas("cGG", "Gamma - Gamma Coin.", size * NCRYSTAL, size * NCRYSTAL); + cGG->Clear();cGG->Divide(NCRYSTAL, NCRYSTAL); - for( int i = 0; i < NCLOVER; i ++){ - for( int j = i+1; j < NCLOVER; j ++){ - cGG->cd( NCLOVER * i + j +1 ); + for( int i = 0; i < NCRYSTAL; i ++){ + for( int j = i+1; j < NCRYSTAL; j ++){ + cGG->cd( NCRYSTAL * i + j +1 ); hgg[i][j]->Draw("colz"); } } @@ -75,10 +104,70 @@ void drawGG(){ } */ +void energyCalibration(int detID = -1, int BG = 10, double threshold = 0.1, double sigmaMax = 5, int peakDensity = 10){ + + TCanvas *cCal = (TCanvas *) gROOT->FindObjectAny("cCal"); + if( cCal == NULL ) cCal = new TCanvas("cCal", "Energy Calibration", 1000, 0, 1000, 600); + cCal->Clear(); -void energyCalibration(){ -// TCanvas + cCal->Divide(2,1); + cCal->SetGrid(); + vector refEnergy = {121.738, + 244.699, + 344.281, + 411.115, + 443.965, + 778.903, + 867.390, + 964.055, + 1085.842, + 1089.700, + 1112.087, + 1408.022}; + double a0[NCRYSTAL]; + double a1[NCRYSTAL]; + + for( int i = 0 ; i < NCRYSTAL; i++){ + if( detID >= 0 && i != detID ) continue; + + cCal->cd(1); + he[i]->Draw(); + vector peaks = fitAuto(he[i], BG, threshold, sigmaMax, peakDensity); + vector> output = FindMatchingPair(peaks, refEnergy); + + vector haha1 = output[0]; + vector haha2 = output[1]; + + TGraph * graph = new TGraph(haha1.size(), &haha1[0], &haha2[0] ); + cCal->cd(2); + graph->Draw("A*"); + + TF1 * fit = new TF1("fit", "pol1" ); + graph->Fit("fit", ""); + + a0[i] = fit->GetParameter(0); + a1[i] = fit->GetParameter(1); + + if( detID < 0 ) { + printf("%2d | a0: %14.10f, a1: %14.10f\n", i, a0[i], a1[i]); + }else{ + printf("%2d | a0, a1 = %14.10f\t%14.10f\n", i, a0[i], a1[i]); + } + } + + if( detID < 0 ){ + FILE * paraOut; + TString filename; + filename.Form("correction_e_auto.dat"); + paraOut = fopen (filename.Data(), "w+"); + printf("=========== save e-correction parameters to %s \n", filename.Data()); + for( int i = 0; i < NCRYSTAL; i++){ + fprintf(paraOut, "%14.10f\t%14.10f\n", a0[i], a1[i]); + } + fflush(paraOut); + fclose(paraOut); + } } diff --git a/correction_e.dat b/correction_e.dat new file mode 100644 index 0000000..3f491ab --- /dev/null +++ b/correction_e.dat @@ -0,0 +1,36 @@ + -0.11851105 0.73602844 + -0.05329699 0.30617141 + 0.21919925 0.30256406 + 0.19549630 0.29988911 + 0.23576185 0.31232252 + 0.15862498 0.31719203 + 0.17456490 0.30495670 + 0.07937138 0.31393499 + -0.17752085 0.30734768 + 0.47543250 0.30864178 + 0.08375230 0.30843852 + 0.28037171 0.31263324 + -0.04410183 0.74143159 + 0.07905962 0.73641543 + -0.05892825 0.71786084 + 0.07476386 0.36785529 + 0.07951184 0.26260823 + 0.02161385 0.25884364 + 0.25371149 0.29681333 + 0.23290589 0.34255969 + 0.20677949 0.30504413 + 0.16341696 0.30761409 + 0.04406586 0.30595347 + 0.07292338 0.30758425 + 0.00136881 0.17925931 + 0.23758200 0.31520725 + -0.47281914 0.17676788 + 0.04230014 0.17917457 + 0.20654489 0.30340233 + 0.20762807 0.30960594 + 0.19673688 0.30110502 + 0.07362825 0.29715807 + 0.17023147 0.30259114 + 0.23642061 0.29846387 + 0.15627111 0.31411607 + -0.01255732 0.15930220 diff --git a/correction_e_auto.dat b/correction_e_auto.dat new file mode 100644 index 0000000..cd64229 --- /dev/null +++ b/correction_e_auto.dat @@ -0,0 +1,36 @@ + -0.11851105 0.73602844 + -0.00614413 0.30614471 + 0.24650404 0.30255455 + 0.11051358 0.29991329 + 0.25582558 0.31231509 + 0.17781452 0.31718544 + 0.16689634 0.30496368 + 0.03831077 0.31393841 + -0.14623990 0.30734342 + 0.35972740 0.30867079 + 0.10284998 0.30843610 + 0.23571282 0.31264887 + -0.04410183 0.74143159 + 0.07905962 0.73641543 + -11.31371425 0.72535320 + 0.06260484 0.36785291 + 0.07951958 0.26260898 + 0.04175235 0.25884441 + 0.16910336 0.29683177 + 0.19876657 0.34256111 + 0.20474761 0.30504873 + 0.07118067 0.30762606 + -0.02019236 0.30596997 + 0.07693261 0.30756913 + 0.00136961 0.17925931 + 0.24935174 0.31520726 + -0.47281914 0.17676788 + 0.02478744 0.17917572 + 0.23969611 0.30339928 + 0.19174182 0.30961199 + 0.17180469 0.30111457 + -0.00984489 0.29718994 + 0.16313776 0.30259590 + 0.20670406 0.29847633 + 0.09021853 0.31413251 + -0.04099465 0.15930879 diff --git a/e_corr.txt b/e_corr.txt deleted file mode 100644 index 1646e28..0000000 --- a/e_corr.txt +++ /dev/null @@ -1,64 +0,0 @@ -0 1 0 0 -0 2 0 0 -1000 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 -0 1 0 0 diff --git a/mapping.h b/mapping.h index 9d63b7c..355d6b7 100644 --- a/mapping.h +++ b/mapping.h @@ -7,7 +7,7 @@ Other : 200 - 299 //==================== mapping -#define NCLOVER 36 +#define NCRYSTAL 36 #define NBGO 9 #define NOTHER 52