Index: /branches/eam_branches/ipp-20110710/psphot/src/psphotPetrosianStats.c
===================================================================
--- /branches/eam_branches/ipp-20110710/psphot/src/psphotPetrosianStats.c	(revision 32315)
+++ /branches/eam_branches/ipp-20110710/psphot/src/psphotPetrosianStats.c	(revision 32316)
@@ -143,15 +143,8 @@
     // if we failed to reach the PETROSIAN_RATIO, use the lowest significant ratio instead (flag this!)
     if (!anyPetro) {
-	// interpolate Rvec between i-1 and i to PETROSIAN_RATIO to get flux (Fvec) and radius (rvec)
-	if (lowestSignificantRadius == 0) { 
-	    // assume Fmax @ R = 0.0
-	    petRadius    = InterpolateValues     (1.0, 0.0, petRatio->data.F32[lowestSignificantRadius], refRadius->data.F32[lowestSignificantRadius], PETROSIAN_RATIO);
-	    petRadiusErr = InterpolateValuesErrX (1.0, 0.0, petRatio->data.F32[lowestSignificantRadius], refRadius->data.F32[lowestSignificantRadius], PETROSIAN_RATIO, 0.0, petRatioErr->data.F32[lowestSignificantRadius]);
-
-	} else {
-	    int n0 = lowestSignificantRadius-1;
-	    int n1 = lowestSignificantRadius;
-	    petRadius    = InterpolateValues     (petRatio->data.F32[n0], refRadius->data.F32[n0], petRatio->data.F32[n1], refRadius->data.F32[n1], PETROSIAN_RATIO);
-	    petRadiusErr = InterpolateValuesErrX (petRatio->data.F32[n0], refRadius->data.F32[n0], petRatio->data.F32[n1], refRadius->data.F32[n1], PETROSIAN_RATIO, petRatioErr->data.F32[n0], petRatioErr->data.F32[n1]);
+	petRadius = refRadius->data.F32[lowestSignificantRadius];
+	petRadiusErr = NAN;
+	if (!isfinite(petRadius)) {
+	    fprintf (stderr, "nan pet radius\n");
 	}
     }
@@ -170,4 +163,7 @@
 		petArea    = InterpolateValues     (refRadius->data.F32[i-1], areaSum->data.F32[i-1], refRadius->data.F32[i], areaSum->data.F32[i], apRadius);
 		petApix    = InterpolateValues     (refRadius->data.F32[i-1], apixSum->data.F32[i-1], refRadius->data.F32[i], apixSum->data.F32[i], apRadius);
+		if (!isfinite(petFlux)) {
+		    fprintf (stderr, "nan pet flux\n");
+		}
 		break;
 	    }
@@ -188,6 +184,7 @@
 	if (!found50 && (fluxSum->data.F32[i] > flux50)) {
 	    if (i == 0) {
-		psWarning ("does this case make any sense? (fluxSum[0] %f > flux50 %f)", fluxSum->data.F32[i], flux50);
-		continue;
+		R50    = InterpolateValues     (fluxSum->data.F32[i], refRadius->data.F32[i], fluxSum->data.F32[i+1], refRadius->data.F32[i+1], flux50);
+		R50err = InterpolateValuesErrX (fluxSum->data.F32[i], refRadius->data.F32[i], fluxSum->data.F32[i+1], refRadius->data.F32[i+1], flux50, sqrt(fluxSumErr2->data.F32[i]), sqrt(fluxSumErr2->data.F32[i+1]));
+		found50 = true;
 	    } else {
 		R50    = InterpolateValues     (fluxSum->data.F32[i-1], refRadius->data.F32[i-1], fluxSum->data.F32[i], refRadius->data.F32[i], flux50);
@@ -198,6 +195,7 @@
 	if (!found90 && (fluxSum->data.F32[i] > flux90)) {
 	    if (i == 0) {
-		psWarning ("does this case make any sense? (fluxSum[0] %f > flux90 %f)", fluxSum->data.F32[i], flux90);
-		continue;
+		R90    = InterpolateValues     (fluxSum->data.F32[i], refRadius->data.F32[i], fluxSum->data.F32[i+1], refRadius->data.F32[i+1], flux90);
+		R90err = InterpolateValuesErrX (fluxSum->data.F32[i], refRadius->data.F32[i], fluxSum->data.F32[i+1], refRadius->data.F32[i+1], flux90, sqrt(fluxSumErr2->data.F32[i]), sqrt(fluxSumErr2->data.F32[i+1]));
+		found90 = true;
 	    } else {
 		R90    = InterpolateValues     (fluxSum->data.F32[i-1], refRadius->data.F32[i-1], fluxSum->data.F32[i], refRadius->data.F32[i], flux90);