After observing the Higgs boson by the ATLAS and CMS experiments at the LHC, accurate measurements of its properties, which allow us to study the electroweak symmetry breaking mechanism, become a high priority for particle physics. The most promising of extracting the Higgs self-coupling at hadron colliders is by examining the double Higgs production, especially in the $b bar{b} gamma gamma$ channel. In this work, we presented full loop calculation for both SM and New Physics effects of the Higgs pair production to next-to-leading-order (NLO), including loop-induced processes $ggto HH$, $ggto HHg$, and $qg to qHH$. We also included the calculation of the corrections from diagrams with only one QCD coupling in $qg to qHH$, which was neglected in the previous studies. With the latest observed limit on the HH production cross-section, we studied the constraints on the effective Higgs couplings for the LHC at center-of-mass energies of 14 TeV and a provisional 100 TeV proton collider within the Future-Circular-Collider (FCC) project. To obtain results better than using total cross-section alone, we focused on the $b bar{b} gamma gamma$ channel and divided the differential cross-section into low and high bins based on the total invariant mass and $p_{T}$ spectra. The new physics effects are further constrained by including extra kinematic information. However, some degeneracy persists, as shown in previous studies, especially in determining the Higgs trilinear coupling. Our analysis shows that the degeneracy is reduced by including the full NLO corrections.