''' @brief Reads in BayesImage (EMC2) log file; returns lists of LogProb,MS,Scale; writes them out, if wished. @author A. Connors ''' # @file get_LogProbMSBScaleKappa_from_log.py # # This script is explicitly for reading in "log" format files # from CHASC EMC2 imaging software; storing values from each iteration # in handy lists for later manipulation; and writing them out in # handy R-Readable txt format. # #from numarray import * ## Changed to numpy Aug 2007 from numpy import * import pyfits ### ### Using the c-format info from DvD's c-code to parse the bayes log-file ### Into something python and/or R can read in as a data-object: ### # #printf("ITERATION NUMBER %d.\n", cont.iter); # 01234567890123456 #printf("Current Log-Posterior: %10g", llike[0].cur); # 01234567890123456789012 #if(cont[0].iter > 1) # printf(" Step Size: %14.10g", llike[0].cur - llike[0].pre); # ### ### right now, ignoring these alphas in results-plot-files: ### printf("Smoothing Hyper-parameters:");printf(" %f ", ms[0].al[level]);printf("\n"); ### # #printf("Expected Total MS Cnt: %g\n", ms[0].ag[ms[0].power][0][0]); # 01234567890123456789012 #printf("Backgound Scale: %g.\n", bkg_scale.scale); # # # #---------OK. SO That was the original C.--------------------------------------# #---------Now I will use it to read and parse the log-file.--------------------# # # # ?There could also be a method to accumulate more into this LogData object??? # # class BayesEMC2LogData: def __init__(self,InLogFile, NumRunIters=1000, NumBurnIn=100, OutDatFile='TestRFormatOutFile.delim'): # # 0/ Initialize handy internal variables: QBkgScale = 'NO' # # 1/ Open log-file in text, read-only format; # WholeLogFileIn = open(InLogFile,'r') # Read in "header" first: try: InputReadString = WholeLogFileIn.readline() except: InputReadString = 'EOF' return # print '**Now At**: ', InputReadString check = '' while(InputReadString!= 'EOF' and check != 'GotToFirstIter'): # 1.1/ Is it bkg scale fit, or not? # (For now, ignores rest of starting info in log) # ie does InputReadString[0:47] = 'A scale parameter will be fit to the bkg model.' # print'Huh! Whats wrong? Length of string: ', len(InputReadString), '\n' if(len(InputReadString) > 47): # print '[0][1][2][3][4][5]etc', InputReadString[0:1],InputReadString[1:2],InputReadString[2:3],InputReadString[3:4],InputReadString[4:5],InputReadString[5:6],InputReadString[6:7],InputReadString[7:8] if (InputReadString[0:47] == 'A scale parameter will be fit to the bkg model.'): QBkgScale = 'YES' print '**Fitting BKGSCALE**' if((len(InputReadString) > 9) and (InputReadString[0:9] != 'ITERATION')): try: InputReadString = WholeLogFileIn.readline() except: InputReadString = 'EOF' break print 'INTRO: ', InputReadString, check = '' elif(len(InputReadString) > 9) and (InputReadString[0:9] == 'ITERATION'): print 'INTRO: ', InputReadString, check = 'GotToFirstIter' else: try: InputReadString = WholeLogFileIn.readline() except: InputReadString = 'EOF' break print 'INTRO: ', InputReadString, # end of check for first 'ITER' # print 'end of check for first ITER' # end of "while", transition to reading first iteration # print 'end of "while", transition to reading first iteration' # # 2/ Get starting iteration; Create/Initialize lists to hold results per iter: # ie InputReadString[0:19] = 'ITERATION NUMBER 1.' if((len(InputReadString) > 19) and (InputReadString[0:18] != 'ITERATION NUMBER 1.')): print '\n Ooops! Some kind of goof! 1st iteration is:', InputReadString,'\n' IterNumber = [1] linesimple = 1 # # # 3/ Then, for each iteration, until the end-of-file: while( (linesimple <= NumRunIters) and (InputReadString != 'EOF')): # print 'Iteration:', linesimple # Then read: Log-Posterior try: InputReadString = WholeLogFileIn.readline() except: InputReadString = 'EOF' break # print '**Reading Log Posterior**: ', InputReadString,InputReadString[23:33] try: InLogPost = float(InputReadString[23:33]) except: print 'Problems with:',InputReadString InputReadString = 'EOF' linesimple = linesimple-1 break # Skip Smoothing hyper-params try: InputReadString = WholeLogFileIn.readline() except: InputReadString = 'EOF' break # print '**Skipping Hyper-Params**: ', InputReadString # Read Expected Total MS Cnt, Backgound Scale try: InputReadString = WholeLogFileIn.readline() except: InputReadString = 'EOF' break # print '**Reading MS line**: ', InputReadString try: InMS = float(InputReadString[23:]) except: print 'Problems with:',InputReadString linesimple = linesimple-1 InputReadString = 'EOF' break if (QBkgScale == 'YES'): try: InputReadString = WholeLogFileIn.readline() except: InputReadString = 'EOF' break # print '**Reading BScale line**: ', InputReadString BScaleString = InputReadString[17:] ## print 'BScale: len, values:', len(BScaleString),',',BScaleString[0:1],',',BScaleString[0:1],',',BScaleString[1:2],',',BScaleString[2:3],',',BScaleString[3:4],',',BScaleString[4:5],',',BScaleString[5:6],',',BScaleString[6:7],',',BScaleString[7:8] try: InBkgScale = float(BScaleString[0:len(BScaleString)-2]) except: print 'Problems with:', BScaleString linesimple = linesimple-1 InputReadString = 'EOF' break # Put them in: if(InputReadString != 'EOF'): if(linesimple == 1): LogPost = [InLogPost] ExpectedMSCounts = [InMS] if (QBkgScale == 'YES'): ExpectedBkgScale = [InBkgScale] # # 4or5/ Open outputfile; write appropriate header # Append these to data-object # And write them to OutFile in appropriate format: elif(linesimple > 1): IterNumber.append(linesimple) LogPost.append(InLogPost) ExpectedMSCounts.append(InMS) if (QBkgScale == 'YES'): ExpectedBkgScale.append(InBkgScale) # Read iteration try: InputReadString = WholeLogFileIn.readline() except: InputReadString = 'EOF' break # print '**Reading ITER line**: ', InputReadString linesimple = linesimple+1 # print 'Incrementing Linesimple:', linesimple, '\n' # # 5/ At end of file, close the output files. # # end of while-not-end-of-file WholeLogFileIn.close() self.IterNumber = IterNumber self.LogPosterior = LogPost self.ExpectedMSCounts = ExpectedMSCounts if (QBkgScale == 'YES'): self.ExpectedBkgScale = ExpectedBkgScale # # 6/ find mean, variance, etc of each variable? # and print it to the screen? # av_LogPost = average(LogPost[NumBurnIn:]) av_ExpectedMSCounts = average(ExpectedMSCounts[NumBurnIn:]) av_ExpectedBkgScale = average(ExpectedBkgScale[NumBurnIn:]) # dl_LogPost = LogPost - av_LogPost # dl_ExpectedMSCounts = ExpectedMSCounts - av_ExpectedMSCounts # dl_ExpectedBkgScale = ExpectedBkgScale -av_ExpectedBkgScale # vr_LogPost = average(dl_LogPost*dl_LogPost) # vr_ExpectedMSCounts = average(dl_ExpectedMSCounts*dl_ExpectedMSCounts) # vr_ExpectedBkgScale = average(dl_ExpectedBkgScale*dl_ExpectedBkgScale) # sqrt_vr_LogPost = sqrt(vr_LogPost) # sqrt_vr_ExpectedMSCounts = sqrt(vr_ExpectedMSCounts) # sqrt_vr_ExpectedBkgScale = sqrt(vr_ExpectedBkgScale) # sk_LogPost = average(dl_LogPost*dl_LogPost)/sqrt_vr_LogPost/sqrt_vr_LogPost/sqrt_vr_LogPost # sk_ExpectedMSCounts = average(dl_ExpectedMSCounts*dl_ExpectedMSCounts*dl_ExpectedMSCounts)/sqrt_vr_ExpectedMSCounts/sqrt_vr_ExpectedMSCounts/sqrt_vr_ExpectedMSCounts # sk_ExpectedBkgScale = average(dl_ExpectedBkgScale*dl_ExpectedBkgScale*dl_ExpectedBkgScale)/sqrt_vr_ExpectedMSCounts/sqrt_vr_ExpectedMSCounts/sqrt_vr_ExpectedMSCounts print '\n Mean, SqrtVar, Skew of LogPosterior : ', av_LogPost #,sqrt_vr_LogPost , sk_LogPost print '\n Mean, SqrtVar, Skew of ExpectedMSCounts: ', av_ExpectedMSCounts#,sqrt_vr_ExpectedMSCounts,sk_ExpectedMSCounts if (QBkgScale == 'YES'): print '\n Mean, SqrtVar, Skew of ExpectedBkgScale: ', av_ExpectedBkgScale#,sqrt_vr_ExpectedBkgScale,sk_ExpectedBkgScale # # 7/ NOW write them out in "Rdelim" format: # testoutfilename = OutDatFile TestRFormatOutFile = open(testoutfilename,'w') # Write out header: delim = ' ' outstring = delim+'LogPostProb'+delim+'ExpectMSCnts' if(QBkgScale == 'YES'): outstring = outstring+delim+'BkgScale' outstring = outstring+'\n' TestRFormatOutFile.write(outstring) for line in range(linesimple-1): outstring = str(IterNumber[line])+delim+str(LogPost[line])+delim+str(ExpectedMSCounts[line]) if(QBkgScale == 'YES'): outstring = outstring+delim+str(ExpectedBkgScale[line]) outstring = outstring+'\n' # print line, outstring TestRFormatOutFile.write(outstring) # end loopover all lines/iterations #------------------------------------------------------------------------------#