|
633 | 633 | "id": "threshold", |
634 | 634 | "metadata": {}, |
635 | 635 | "outputs": [], |
636 | | - "source": [ |
637 | | - "def estimate_threshold_crossings(results, distances):\n", |
638 | | - " \"\"\"Estimate threshold from crossings between consecutive distance curves.\"\"\"\n", |
639 | | - " crossings = []\n", |
640 | | - "\n", |
641 | | - " for idx in range(len(distances) - 1):\n", |
642 | | - " d1 = distances[idx]\n", |
643 | | - " d2 = distances[idx + 1]\n", |
644 | | - "\n", |
645 | | - " pts1 = {p: ler for (d, p, ler, _) in results if d == d1}\n", |
646 | | - " pts2 = {p: ler for (d, p, ler, _) in results if d == d2}\n", |
647 | | - "\n", |
648 | | - " common_ps = sorted(set(pts1.keys()) & set(pts2.keys()))\n", |
649 | | - " if len(common_ps) < 2:\n", |
650 | | - " continue\n", |
651 | | - "\n", |
652 | | - " for j in range(len(common_ps) - 1):\n", |
653 | | - " p_lo = common_ps[j]\n", |
654 | | - " p_hi = common_ps[j + 1]\n", |
655 | | - " diff_lo = pts1[p_lo] - pts2[p_lo]\n", |
656 | | - " diff_hi = pts1[p_hi] - pts2[p_hi]\n", |
657 | | - "\n", |
658 | | - " if diff_lo * diff_hi < 0:\n", |
659 | | - " t = diff_lo / (diff_lo - diff_hi)\n", |
660 | | - " p_cross = p_lo + t * (p_hi - p_lo)\n", |
661 | | - " crossings.append((d1, d2, p_cross))\n", |
662 | | - "\n", |
663 | | - " return crossings\n", |
664 | | - "\n", |
665 | | - "\n", |
666 | | - "def estimate_threshold_fit(results, distances, fit_func=func, p0=None):\n", |
667 | | - " \"\"\"Estimate threshold using universal scaling fit.\"\"\"\n", |
668 | | - " filtered = [(d, p, ler) for (d, p, ler, _) in results if ler > 0 and d in distances]\n", |
669 | | - " if len(filtered) < 5:\n", |
670 | | - " return None\n", |
671 | | - "\n", |
672 | | - " plist = np.array([p for _, p, _ in filtered])\n", |
673 | | - " dlist = np.array([float(d) for d, _, _ in filtered])\n", |
674 | | - " plog = np.array([ler for _, _, ler in filtered])\n", |
675 | | - "\n", |
676 | | - " if p0 is None:\n", |
677 | | - " if fit_func == func6:\n", |
678 | | - " p0 = [0.5, 0.005]\n", |
679 | | - " else:\n", |
680 | | - " p0 = [0.005, 1.5, 0.5, 1.0, 0.0]\n", |
681 | | - "\n", |
682 | | - " try:\n", |
683 | | - " popt, stdev = threshold_fit(plist, dlist, plog, fit_func, p0)\n", |
684 | | - " return popt, stdev\n", |
685 | | - " except Exception as e:\n", |
686 | | - " print(f\" Fit failed: {e}\")\n", |
687 | | - " return None\n", |
688 | | - "\n", |
689 | | - "\n", |
690 | | - "# === Build comparison table ===\n", |
691 | | - "print(\"=== Threshold Comparison ===\")\n", |
692 | | - "print()\n", |
693 | | - "header = f\"{'Variant':<22} | {'Z crossing':>10} | {'Z fit (func6)':>16} | {'X crossing':>10} | {'X fit (func6)':>16}\"\n", |
694 | | - "print(header)\n", |
695 | | - "print(\"-\" * len(header))\n", |
696 | | - "\n", |
697 | | - "for v in VARIANTS:\n", |
698 | | - " row = f\"{variant_labels[v]:<22}\"\n", |
699 | | - "\n", |
700 | | - " for basis_results in [results_z, results_x]:\n", |
701 | | - " data = basis_results[v]\n", |
702 | | - "\n", |
703 | | - " # Crossing estimate\n", |
704 | | - " crossings = estimate_threshold_crossings(data, DISTANCES)\n", |
705 | | - " if crossings:\n", |
706 | | - " avg_cross = np.mean([c[2] for c in crossings])\n", |
707 | | - " row += f\" | {avg_cross:>10.4f}\"\n", |
708 | | - " else:\n", |
709 | | - " row += f\" | {'N/A':>10}\"\n", |
710 | | - "\n", |
711 | | - " # Fit estimate\n", |
712 | | - " result = estimate_threshold_fit(data, DISTANCES, func6, p0=[0.5, 0.005])\n", |
713 | | - " if result is not None:\n", |
714 | | - " popt, stdev = result\n", |
715 | | - " row += f\" | {popt[1]:>7.4f}+/-{stdev[1]:<5.4f}\"\n", |
716 | | - " else:\n", |
717 | | - " row += f\" | {'N/A':>16}\"\n", |
718 | | - "\n", |
719 | | - " print(row)\n", |
720 | | - "\n", |
721 | | - "print()\n", |
722 | | - "\n", |
723 | | - "# === Detailed crossing info ===\n", |
724 | | - "print(\"\\n=== Detailed Curve Crossings ===\")\n", |
725 | | - "for v in VARIANTS:\n", |
726 | | - " print(f\"\\n--- {variant_labels[v]} ---\")\n", |
727 | | - " for basis, basis_results in [(\"Z\", results_z[v]), (\"X\", results_x[v])]:\n", |
728 | | - " crossings = estimate_threshold_crossings(basis_results, DISTANCES)\n", |
729 | | - " print(f\" {basis} basis:\")\n", |
730 | | - " if crossings:\n", |
731 | | - " for d1, d2, p_cross in crossings:\n", |
732 | | - " print(f\" d={d1}/d={d2} crossing at p ~ {p_cross:.4f}\")\n", |
733 | | - " avg = np.mean([c[2] for c in crossings])\n", |
734 | | - " print(f\" Average: p ~ {avg:.4f}\")\n", |
735 | | - " else:\n", |
736 | | - " print(\" No crossings found.\")" |
737 | | - ] |
| 636 | + "source": "def estimate_threshold_crossings(results, distances):\n \"\"\"Estimate threshold from crossings between consecutive distance curves.\"\"\"\n crossings = []\n\n for idx in range(len(distances) - 1):\n d1 = distances[idx]\n d2 = distances[idx + 1]\n\n pts1 = {p: ler for (d, p, ler, _) in results if d == d1}\n pts2 = {p: ler for (d, p, ler, _) in results if d == d2}\n\n common_ps = sorted(set(pts1.keys()) & set(pts2.keys()))\n if len(common_ps) < 2:\n continue\n\n for j in range(len(common_ps) - 1):\n p_lo = common_ps[j]\n p_hi = common_ps[j + 1]\n diff_lo = pts1[p_lo] - pts2[p_lo]\n diff_hi = pts1[p_hi] - pts2[p_hi]\n\n if diff_lo * diff_hi < 0:\n t = diff_lo / (diff_lo - diff_hi)\n p_cross = p_lo + t * (p_hi - p_lo)\n crossings.append((d1, d2, p_cross))\n\n return crossings\n\n\ndef estimate_threshold_fit(results, distances, fit_func=func, p0=None):\n \"\"\"Estimate threshold using universal scaling fit.\"\"\"\n filtered = [(d, p, ler) for (d, p, ler, _) in results if ler > 0 and d in distances]\n if len(filtered) < 5:\n return None\n\n plist = np.array([p for _, p, _ in filtered])\n dlist = np.array([float(d) for d, _, _ in filtered])\n plog = np.array([ler for _, _, ler in filtered])\n\n if p0 is None:\n p0 = [0.5, 0.005] if fit_func == func6 else [0.005, 1.5, 0.5, 1.0, 0.0]\n\n try:\n popt, stdev = threshold_fit(plist, dlist, plog, fit_func, p0)\n except (RuntimeError, ValueError) as e:\n print(f\" Fit failed: {e}\")\n return None\n else:\n return popt, stdev\n\n\n# === Build comparison table ===\nprint(\"=== Threshold Comparison ===\")\nprint()\nheader = f\"{'Variant':<22} | {'Z crossing':>10} | {'Z fit (func6)':>16} | {'X crossing':>10} | {'X fit (func6)':>16}\"\nprint(header)\nprint(\"-\" * len(header))\n\nfor v in VARIANTS:\n row = f\"{variant_labels[v]:<22}\"\n\n for basis_results in [results_z, results_x]:\n data = basis_results[v]\n\n # Crossing estimate\n crossings = estimate_threshold_crossings(data, DISTANCES)\n if crossings:\n avg_cross = np.mean([c[2] for c in crossings])\n row += f\" | {avg_cross:>10.4f}\"\n else:\n row += f\" | {'N/A':>10}\"\n\n # Fit estimate\n result = estimate_threshold_fit(data, DISTANCES, func6, p0=[0.5, 0.005])\n if result is not None:\n popt, stdev = result\n row += f\" | {popt[1]:>7.4f}+/-{stdev[1]:<5.4f}\"\n else:\n row += f\" | {'N/A':>16}\"\n\n print(row)\n\nprint()\n\n# === Detailed crossing info ===\nprint(\"\\n=== Detailed Curve Crossings ===\")\nfor v in VARIANTS:\n print(f\"\\n--- {variant_labels[v]} ---\")\n for basis, basis_results in [(\"Z\", results_z[v]), (\"X\", results_x[v])]:\n crossings = estimate_threshold_crossings(basis_results, DISTANCES)\n print(f\" {basis} basis:\")\n if crossings:\n for d1, d2, p_cross in crossings:\n print(f\" d={d1}/d={d2} crossing at p ~ {p_cross:.4f}\")\n avg = np.mean([c[2] for c in crossings])\n print(f\" Average: p ~ {avg:.4f}\")\n else:\n print(\" No crossings found.\")" |
738 | 637 | }, |
739 | 638 | { |
740 | 639 | "cell_type": "markdown", |
|
0 commit comments