CHICOS Maximum Likelihood Code

CTShower Code (May not be most recent.)
Shower Test Page (V. 3.03 Testing may not be most current. Be wary of using v. 3.02 as you may have to kill it manually.)
Write-up Updated.

From looking at dummy reconstructions, it seems that there is something wrong in the way the intensity probabilities are being computed. Higher energies are giving a lower likelihood even when the LDF appears to match better (and the chisquare is better).

Looking at the actual numbers, something does appear to be wrong. I had the program print out the numbers going into the calculation P_I = P_Poisson(Total I; LDF). Two things are immediately apparent: there is only one calculation per site (there should be one per shmoo) and the LDF numbers do not match up with the "fit" values when one does a dummy reconstruction on the web with the same energy.

02-27-07, E=1.0e19
Energy = 1e+19
Total I = 3.91827 LDF = 4.97659 PI = -1.74447
Total I = 0.908654 LDF = 0.225396 PI = -1.54332
Total I = 2.98361 LDF = 3.48503 PI = -1.53127
Total I = 1.21359 LDF = 0.739656 PI = -1.21004
Total I = 0.322581 LDF = 0.0289634 PI = -1.05975
02-27-07, E=1.5e19
nergy = 1.5e+19
Total I = 3.91827 LDF = 7.34132 PI = -2.58588
Total I = 0.908654 LDF = 0.332967 PI = -1.29635
Total I = 2.98361 LDF = 5.14037 PI = -2.02704
Total I = 1.21359 LDF = 1.09224 PI = -1.08956
Total I = 0.322581 LDF = 0.0427982 PI = -0.947625

The numbers that appear when you do a dummy reconstruction online are the measured intensity in #/m^2 (averaged over both shmoos at a site) and the value of the LDF (in #/m^2) at that location.

Going by the code, the numbers that should be going into the P_I calculation are the total number of particles collected by each shmoo, and the LDF*Area/cos(theta). The actual "totalInten" however, appears to simply be the #/m^2 averaged over both shmoos. I'm not sure what the "expected_inten" is.

Problems:
- Total intensity should be shmoo-by-shmoo, not site-by-site!
- Are we comparing particle count or intensity per meter squared?
- Does hit->GetInten() retrieve particle count or intensity per m^2? Does this vary with version? (See some earlier notes on this.)

Actually, the first problem may be limited to Dummy Reconstructions, which did not use Shawn's fitshower object... Fixing...

OK, now I get the following:

02-27-07, E=1.0e19
Energy = 1e+19
Total I = 5.3 LDF = 4.97659 PI = -1.77878
Total I = 2.85 LDF = 4.97659 PI = -2.00965
Total I = 1.89 LDF = 0.229773 PI = -3.60339
Total I = 3.64 LDF = 3.48503 PI = -1.59116
Total I = 2.5 LDF = 0.739656 PI = -2.69455
Total I = 0.4 LDF = 0.0289634 PI = -1.32604
02-27-07, E=1.5e19
Energy = 1.5e+19
Total I = 5.3 LDF = 7.34132 PI = -2.08301
Total I = 2.85 LDF = 7.34132 PI = -3.26637
Total I = 1.89 LDF = 0.339432 PI = -2.97561
Total I = 3.64 LDF = 5.14037 PI = -1.83182
Total I = 2.5 LDF = 1.09224 PI = -2.07263
Total I = 0.4 LDF = 0.0427982 PI = -1.18369

The likelihood changes slightly, and individual shmoos are now being looked at. It is now clear that the total_Inten is the raw number of hits at a shmoo. The expected_inten should be calculated as the LDF*shmoo_area/cos(theta). The expected_inten does not seem to correspond with what appears as the "fit" value in the web reconstruction output. A problem with the angle?

Testing removal of the cos(theta) denominator...

02-27-07, E=1.0e19
Energy = 1e+19
Total I = 5.3 LDF = 4.97587 PI = -1.77883
Total I = 2.85 LDF = 4.97587 PI = -2.00934
Total I = 1.89 LDF = 0.22974 PI = -3.60363
Total I = 3.64 LDF = 3.48452 PI = -1.59119
Total I = 2.5 LDF = 0.739549 PI = -2.69481
Total I = 0.4 LDF = 0.0289592 PI = -1.32609
02-27-07, E=1.5e19
Energy = 1.5e+19
Total I = 5.3 LDF = 7.34025 PI = -2.08272
Total I = 2.85 LDF = 7.34025 PI = -3.26572
Total I = 1.89 LDF = 0.339383 PI = -2.97583
Total I = 3.64 LDF = 5.13962 PI = -1.8316
Total I = 2.5 LDF = 1.09208 PI = -2.07284
Total I = 0.4 LDF = 0.0427919 PI = -1.18374

Interestingly, the expected_inten did not change all that much. The angle for this shower is high, this should have changed the expected intensity by 40%. What was cos(theta)_up to begin with?

Update - there is a problem with theta. Theta for the shower is close to 0, cos(theta) is close to 1, even though the shower is inclined to 56 degrees... Still need to look into this.