DaubechiesGenerator.java

package com.morphiqlabs.wavelet.util;

/**
 * Utility class for computing Daubechies wavelet coefficients using spectral factorization.
 * 
 * This class implements the classic algorithm for generating Daubechies scaling function
 * coefficients based on the construction from "Ten Lectures on Wavelets" by I. Daubechies.
 * 
 * Supports orders 1-38 (DB2-DB38 in standard notation). DB38 is the maximum order
 * as higher orders suffer from numerical instability issues.
 */
public final class DaubechiesGenerator {

    /**
     * Generate Daubechies scaling coefficients for a given order using the
     * spectral factorization method.
     * 
     * @param N the number of vanishing moments (order)
     * @return the scaling function coefficients h[k]
     */
    public static double[] generateCoefficients(int N) {
        if (N < 1 || N > 38) {
            throw new IllegalArgumentException("Order must be between 1 and 38 (PyWavelets maximum)");
        }
        
        // Get raw coefficients and normalize them
        double[] raw = getRawCoefficients(N);
        return normalizeCoefficients(raw);
    }

    /**
     * Get raw (possibly unnormalized) coefficients.
     */
    private static double[] getRawCoefficients(int N) {
        return switch (N) {
            case 12 -> {
                // DB12 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -0.00000152907175807, 0.00001277695221938, -0.00002424154575703, -0.00008850410920820,
                    0.00038865306282093, 0.00000654512821251, -0.00217950361862776, 0.00224860724099524,
                    0.00671149900879551, -0.01284082519830068, -0.01221864906974828, 0.04154627749508444,
                    0.01084913025582218, -0.09643212009650708, 0.00535956967435215, 0.18247860592757967,
                    -0.02377925725606973, -0.31617845375278553, -0.04476388565377463, 0.51588647842781565,
                    0.65719872257930712, 0.37735513521421266, 0.10956627282118515, 0.01311225795722952
                };
            }
            case 14 -> {
                // DB14 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -0.00000017871399683, 0.00000172499467537, -0.00000438970490178, -0.00001033720918457,
                    0.00006875504252698, -0.00004177724577037, -0.00038683194731295, 0.00070802115423553,
                    0.00106169108560676, -0.00384963886802219, -0.00074621898926838, 0.01278949326633341,
                    -0.00561504953035696, -0.03018535154039063, 0.02698140830791292, 0.05523712625921604,
                    -0.07154895550404614, -0.08674841156816969, 0.13998901658446070, 0.13839521386480660,
                    -0.21803352999327605, -0.27168855227874805, 0.21867068775890652, 0.63118784910485681,
                    0.55430561794089384, 0.25485026779262138, 0.06236475884939890, 0.00646115346008795
                };
            }
            case 16 -> {
                // DB16 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -0.00000002109339630, 0.00000023087840869, -0.00000073636567855, -0.00000104357134231,
                    0.00001133660866128, -0.00001394566898821, -0.00006103596621411, 0.00017478724522534,
                    0.00011424152003872, -0.00094102174935957, 0.00040789698084971, 0.00312802338120627,
                    -0.00364427962149839, -0.00699001456341392, 0.01399376885982873, 0.01029765964095597,
                    -0.03688839769173014, -0.00758897436885774, 0.07592423604427631, -0.00623972275247487,
                    -0.13238830556381040, 0.02734026375271604, 0.21119069394710430, -0.02791820813302828,
                    -0.32706331052791771, -0.08975108940248964, 0.44029025688635692, 0.63735633208378895,
                    0.43031272284600380, 0.16506428348885313, 0.03490771432367334, 0.00318922092534774
                };
            }
            case 18 -> {
                // DB18 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -0.00000000250793445, 0.00000003068835863, -0.00000011760987670, -0.00000007691632690,
                    0.00000176871298363, -0.00000333263447889, -0.00000852060253745, 0.00003741237880740,
                    -0.00000015359171235, -0.00019864855231175, 0.00021358156191034, 0.00062846568296515,
                    -0.00134059629833611, -0.00111873266699250, 0.00494334360546674, 0.00011863003385812,
                    -0.01305148094661200, 0.00626216795430571, 0.02667070592647059, -0.02373321039586000,
                    -0.04452614190298233, 0.05705124773853688, 0.06488721621190545, -0.10675224665982849,
                    -0.09233188415084628, 0.16708131276325741, 0.14953397556537779, -0.21648093400514298,
                    -0.29365404073655876, 0.14722311196992816, 0.57180165488865131, 0.57182680776660721,
                    0.31467894133703173, 0.10358846582242359, 0.01928853172414638, 0.00157631021844076
                };
            }
            case 20 -> {
                // DB20 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -0.00000000029988365, 0.00000000405612706, -0.00000001814843248, 0.00000000020143220,
                    0.00000026339242263, -0.00000068470795970, -0.00000101199401002, 0.00000724124828767,
                    -0.00000437614386218, -0.00003710586183395, 0.00006774280828378, 0.00010153288973670,
                    -0.00038510474869922, -0.00005349759843998, 0.00139255961932314, -0.00083156217282256,
                    -0.00358149425960962, 0.00442054238704579, 0.00672162730225946, -0.01381052613715192,
                    -0.00878932492390156, 0.03229429953076958, 0.00587468181181183, -0.06172289962468046,
                    0.00563224685730744, 0.10229171917444256, -0.02471682733861359, -0.15545875070726795,
                    0.03985024645777120, 0.22829105081991632, -0.01672708830907701, -0.32678680043403496,
                    -0.13921208801148388, 0.36150229873933104, 0.61049323893859386, 0.47269618531090168,
                    0.21994211355139703, 0.06342378045908152, 0.01054939462495040, 0.00077995361366685
                };
            }
            case 22 -> {
                // DB22 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -0.00000000003602113, 0.00000000053359388, -0.00000000272962315, 0.00000000168017140,
                    0.00000003761228749, -0.00000012833362288, -0.00000008779879873, 0.00000129518205732,
                    -0.00000156517913200, -0.00000616672931647, 0.00001737375695756, 0.00001137434966213,
                    -0.00009405223634816, 0.00004345899904532, 0.00032860941421368, -0.00042378739983918,
                    -0.00077069098812312, 0.00182701049565728, 0.00104426073918603, -0.00545569198615672,
                    0.00030013739850764, 0.01256472521834337, -0.00621378284936466, -0.02348000134449319,
                    0.02058670762756536, 0.03697084662069802, -0.04653081182750671, -0.05136425429744413,
                    0.08455737636682607, 0.06807631439273222, -0.13176813768668341, -0.09711079840911471,
                    0.17997318799289130, 0.16409318810676649, -0.20056840610488710, -0.31272658042829621,
                    0.07372450118363015, 0.50790109062216393, 0.57843273100952441, 0.36772868344603749,
                    0.14836754089011142, 0.03806993723641108, 0.00572185463133454, 0.00038626323149110
                };
            }
            case 24 -> {
                // DB24 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -0.00000000000434278, 0.00000000006991801, -0.00000000040246586, 0.00000000047483758,
                    0.00000000515777679, -0.00000002255740388, -0.00000000050576454, 0.00000021663396533,
                    -0.00000040325077569, -0.00000089802531439, 0.00000390110033860, 0.00000001341157751,
                    -0.00002022888292613, 0.00002183241460467, 0.00006559388639306, -0.00014600798177626,
                    -0.00011812332379696, 0.00058612705931831, -0.00004416184856142, -0.00169645681897482,
                    0.00115376493683948, 0.00373604617828252, -0.00474656878632311, -0.00629143537001819,
                    0.01304997087108574, 0.00766172188164659, -0.02821310709490189, -0.00494470942812563,
                    0.05130162003998088, -0.00457843624181922, -0.08216165420800167, 0.02098011370914481,
                    0.12101630346922423, -0.03877717357792002, -0.17117535137034690, 0.04252872964148383,
                    0.23923738878031087, 0.00477661368434473, -0.31794307899936275, -0.18727140688515623,
                    0.28098555323371188, 0.57493922109554196, 0.50437104083992501, 0.27290891606772633,
                    0.09726223583362520, 0.02248233994971641, 0.00308208171490549, 0.00019143580094755
                };
            }
            case 26 -> {
                // DB26 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -5.25187122424443475E-13, 9.13051001637179662E-12, -5.84040818534117088E-11, 1.00230319104652690E-10, 
                    6.78004724582863666E-10, -3.77601047853232407E-09, 2.16932825985032301E-09, 3.40779562129072982E-08, 
                    -8.90446637016859012E-08, -1.07900423757867145E-07, 7.93921063370995239E-07, -4.65046322064026265E-07, 
                    -3.88740016185679531E-06, 7.00007868296498702E-06, 1.07422154087219501E-05, -4.10967399639147754E-05, 
                    -5.27779549303786932E-06, 1.57479523860749350E-04, -1.06057474828380395E-04, -4.31955707426180766E-04, 
                    6.16138220457434441E-04, 8.38348805654361628E-04, -2.14553028156762087E-03, -9.39058250473828954E-04, 
                    5.60194723942380492E-03, -5.28738399262681461E-04, -1.17854979061930293E-02, 5.82958055531888810E-03, 
                    2.07349201799638255E-02, -1.77609035683581849E-02, -3.13781103630677571E-02, 3.85357159711118627E-02, 
                    4.22321857963720362E-02, -6.86547596040359143E-02, -5.34485616814831949E-02, 1.06482405249808634E-01, 
                    6.98231861132923709E-02, -1.47977193275254493E-01, -1.04323900285927043E-01, 1.82755409589672374E-01, 
                    1.81291832311122697E-01, -1.74839961289392498E-01, -3.26384593691780023E-01, 1.77407678098668578E-03, 
                    4.39158311789166256E-01, 5.73669043034222281E-01, 4.13292962278356379E-01, 1.95039438716770097E-01, 
                    6.22747440251496046E-02, 1.30975542925585008E-02, 1.65052023353298824E-03, 9.49379575071059272E-05
                };
            }
            case 28 -> {
                // DB28 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -6.36777235471485717E-14, 1.18885053340590150E-12, -8.36549047125880063E-12, 1.86736726378339058E-11, 
                    8.49222001105638229E-11, -6.07704124722901063E-10, 6.94454032894622681E-10, 5.04404705638343681E-09, 
                    -1.78413869087570998E-08, -8.26238731562655765E-09, 1.49066001353536223E-07, -1.75746117320984265E-07, 
                    -6.67021547995489289E-07, 1.84036373451776915E-06, 1.24790031757483418E-06, -1.00432604133342260E-05, 
                    4.63866498139429476E-06, 3.64140121105080249E-05, -4.90771341619025055E-05, -8.90390149004448772E-05, 
                    2.29579098223345627E-04, 1.15465606365892129E-04, -7.48674955911463004E-04, 1.41567239314046439E-04, 
                    1.87599866820279564E-03, -1.36037384563969239E-03, -3.72546124707425494E-03, 4.78486311245424154E-03, 
                    5.83881662774894498E-03, -1.20635919682184900E-02, -6.81554976455230922E-03, 2.46880600101518667E-02, 
                    4.43173291006298837E-03, -4.33333686160862833E-02, 3.44801895554095123E-03, 6.77478955019093360E-02, 
                    -1.73419228313058983E-02, -9.76853558056524351E-02, 3.44786312750997026E-02, 1.34627567910226092E-01, 
                    -4.68382337445516772E-02, -1.82877330732984927E-01, 3.69068853157112700E-02, 2.45808151373759554E-01, 
                    3.28578791633871020E-02, -3.01327809532641766E-01, -2.30498954047582527E-01, 2.00176144045984439E-01, 
                    5.30516293441485765E-01, 5.24998231630335543E-01, 3.22563361285522432E-01, 1.35137914253641050E-01, 
                    3.90926081154053459E-02, 7.54265037764685880E-03, 8.79498515984386987E-04, 4.71080777501405105E-05
                };
            }
            case 30 -> {
                // DB30 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -7.73794263095440493E-15, 1.54399757084761997E-13, -1.18523759210158225E-12, 3.23942863853228594E-12, 
                    1.00010513139317123E-11, -9.46138799727680259E-11, 1.61362297827090424E-10, 6.98486269183218252E-10, 
                    -3.33110568046757819E-09, 5.55339786139705414E-10, 2.60544275497762536E-08, -4.76437996513945334E-08, 
                    -1.00041468235450085E-07, 4.26166232601157230E-07, 1.09947433852620331E-08, -2.18726767699616650E-06, 
                    2.32754909849368655E-06, 7.25214553589046890E-06, -1.63615247872542659E-05, -1.33971686329397169E-05, 
                    6.98200837080832767E-05, -8.54830546758407029E-06, -2.16171830116963369E-04, 1.72482584235170963E-04, 
                    5.05094823903346792E-04, -7.67878250438091864E-04, -8.60927696811042398E-04, 2.32452009406009918E-03, 
                    8.43384586662093445E-04, -5.53073014819200313E-03, 6.19671756497724405E-04, 1.09156316583048901E-02, 
                    -5.29685966613108704E-03, -1.83997438681173416E-02, 1.52879607698573963E-02, 2.70786195952941837E-02, 
                    -3.22637589193522090E-02, -3.56733974967596082E-02, 5.67123657447356974E-02, 4.38016646714177310E-02, 
                    -8.76586900363836574E-02, -5.38064654582570759E-02, 1.22747746045009376E-01, 7.27786589703644238E-02, 
                    -1.57236817959993808E-01, -1.14558219432707775E-01, 1.77829873244836734E-01, 1.99462121580664314E-01, 
                    -1.41968513330082924E-01, -3.32966975020855593E-01, -6.61836707759373144E-02, 3.66242683371627964E-01, 
                    5.57572232912836419E-01, 4.50487821853317816E-01, 2.42020670940214094E-01, 9.12383040670157047E-02, 
                    2.41308326715883800E-02, 4.30079716504806926E-03, 4.66637950428550909E-04, 2.33861617273142149E-05
                };
            }
            case 32 -> {
                // DB32 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -9.42101913953507886E-16, 2.00071530381052503E-14, -1.66380048943340227E-13, 5.36148222961180208E-13, 
                    1.07561065350106221E-12, -1.43091876516920236E-11, 3.26327074133290788E-11, 8.90472379622160575E-11, 
                    -5.88109146263460551E-10, 4.38438779994047434E-10, 4.25042231198059262E-09, -1.10438302172264895E-08, 
                    -1.21992435948337315E-08, 8.96596631195772875E-08, -5.00336186874823014E-08, -4.28597069315145724E-07, 
                    7.56004762559594814E-07, 1.20288903632162094E-06, -4.55830957626442341E-06, -6.36178153226025547E-07, 
                    1.82426840198069137E-05, -1.29404577940551275E-05, -5.25980928268432311E-05, 8.10367832913483832E-05, 
                    1.05391546173982812E-04, -3.05965442382691191E-04, -1.02453731060739617E-04, 8.67305851845055504E-04, 
                    -2.21167872957909794E-04, -1.96474055582177834E-03, 1.46895510046846783E-03, 3.62722464068786481E-03, 
                    -4.64921675118441183E-03, -5.41156825727579123E-03, 1.10174007154068814E-02, 6.16752731068567501E-03, 
                    -2.16628228363911941E-02, -4.14590766082721836E-03, 3.70514579235446812E-02, -2.38026446493257377E-03, 
                    -5.69263140624784378E-02, 1.41061515161066079E-02, 8.08741406384839573E-02, -2.96278725084477036E-02, 
                    -1.09456113116089382E-01, 4.44049081999397383E-02, 1.45232079475286657E-01, -4.89951171846717409E-02, 
                    -1.92102344708546896E-01, 2.46624448396974040E-02, 2.48310642356880162E-01, 6.47133548055162378E-02, 
                    -2.77421581558427222E-01, -2.66698181476675567E-01, 1.20630538265617829E-01, 4.77809163733948383E-01, 
                    5.34317919340953851E-01, 3.67509628597349647E-01, 1.75750783639438912E-01, 6.02574991203353727E-02, 
                    1.46810463814191355E-02, 2.43126191957226609E-03, 2.46656690638090328E-04, 1.16146330213501490E-05
                };
            }
            case 34 -> {
                // DB34 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -1.14894475448059004E-16, 2.58733838193569964E-15, -2.31708370390640842E-14, 8.57919405179973320E-14, 
                    9.79945115821159831E-14, -2.10787910891530174E-12, 6.08012535400016750E-12, 1.00420873546176984E-11, 
                    -9.90477453763240945E-11, 1.30041031860941529E-10, 6.44637821032340189E-10, -2.31650194699548295E-09, 
                    -8.66574426136872175E-10, 1.74042333293606813E-08, -1.99034650153173689E-08, -7.52670174041258947E-08, 
                    2.02599066666785925E-07, 1.44819570833318502E-07, -1.11630653481700838E-06, 4.97971810142130811E-07, 
                    4.16987175854702820E-06, -5.71082651099830398E-06, -1.05765749425795057E-05, 2.84495141969780752E-05, 
                    1.35311722724964963E-05, -9.91469777078013458E-05, 2.66005001845344188E-05, 2.65077239755805768E-04, 
                    -2.32673214023353158E-04, -5.52735576214419765E-04, 8.75199906407868908E-04, 8.58995987436366164E-04, 
                    -2.39945394353705595E-03, -7.69212797506783723E-04, 5.33495076875993568E-03, -6.19474884515387324E-04, 
                    -1.00455067083615197E-02, 4.71364926099980975E-03, 1.64093741998651912E-02, -1.31439800166571613E-02, 
                    -2.36717379228263657E-02, 2.72283507563541964E-02, 3.07397465739593437E-02, -4.74385596452777600E-02, 
                    -3.70128384178624523E-02, 7.31852354367956009E-02, 4.35760946496312959E-02, -1.02947596992814083E-01, 
                    -5.44829680641390479E-02, 1.34125960271136130E-01, 7.79918469379481116E-02, -1.60924927177866800E-01, 
                    -1.27337358223801156E-01, 1.66601750412207456E-01, 2.16907220187427585E-01, -1.03891915515640476E-01, 
                    -3.31525301508386938E-01, -1.28246842174437159E-01, 2.90366329507274978E-01, 5.30555099656463192E-01, 
                    4.78478746279371037E-01, 2.87765059233714537E-01, 1.24152482111376805E-01, 3.90488413517859415E-02, 
                    8.81988940388497844E-03, 1.36406139005904991E-03, 1.29947620067953005E-04, 5.77051063273028522E-06
                };
            }
            case 36 -> {
                // DB36 coefficients from PyWavelets (already normalized)
                yield new double[]{
                    -1.40327417537319069E-17, 3.33997198481869300E-16, -3.20462854340174970E-15, 1.33807138629910590E-14, 
                    5.54226318263980437E-15, -3.02928502697487722E-13, 1.07096935711401707E-12, 8.87684628721737460E-13, 
                    -1.59971668926135704E-11, 3.03742909811253504E-11, 8.96241820385961248E-11, -4.51254577856324943E-10, 
                    1.09081555371375182E-10, 3.13884169578242403E-09, -5.61278434332779106E-09, -1.15609368881700846E-08, 
                    4.79904346545099230E-08, 2.75324907333951215E-09, -2.45537765843423269E-07, 2.54842352255657764E-07, 
                    8.31142127970777899E-07, -1.87081160285918077E-06, -1.58614578243457745E-06, 8.37221819816078826E-06, 
                    -1.18347105998561593E-06, -2.73139082465433775E-05, 2.37510668366086083E-05, 6.69474119693058993E-05, 
                    -1.13189946808466566E-04, -1.15511889584352714E-04, 3.69350728496751048E-04, 8.61456575899270194E-05, 
                    -9.46340382326110166E-04, 2.77681279571202614E-04, 1.99079377185173729E-03, -1.50307406629664382E-03, 
                    -3.48454144540488340E-03, 4.41348483535057567E-03, 5.02298910666582922E-03, -9.99026347328137165E-03, 
                    -5.65781324505881811E-03, 1.90635947806253592E-02, 3.98404019871700498E-03, -3.19807206776396985E-02, 
                    1.42497266176539166E-03, 4.85130835478090883E-02, -1.13191003168174285E-02, -6.82090166368175127E-02, 
                    2.50387214495684900E-02, 9.11567822580165443E-02, -3.98808535755131727E-02, -1.18803754310135637E-01, 
                    5.02761800735384290E-02, 1.54106236627642890E-01, -4.58614007463927151E-02, -1.99337205608649620E-01, 
                    7.27851509579222934E-03, 2.46537277608974204E-01, 9.81142041631147682E-02, -2.46807036978125532E-01, 
                    -2.94421039589114586E-01, 4.39751975293486280E-02, 4.17875335600969788E-01, 5.32266895260728679E-01, 
                    4.06433697708255326E-01, 2.17756953097900802E-01, 8.56520925952640871E-02, 2.48905656448279652E-02, 
                    5.24029737740988426E-03, 7.60215109966848844E-04, 6.82602867854635820E-05, 2.86792518275594617E-06
                };
            }
            case 38 -> {
                // DB38 coefficients from PyWavelets (already normalized) - PyWavelets maximum
                yield new double[]{
                    -1.71615245108874421E-18, 4.30459683955879032E-17, -4.40530704248346116E-16, 2.04509967678898872E-15, 
                    -4.56339716212737354E-16, -4.24981781957146316E-14, 1.80866123627453061E-13, 2.62649650406525195E-14, 
                    -2.48478923756364274E-12, 6.29153731703950837E-12, 1.10169293459945451E-11, -8.27825652253813440E-11, 
                    6.73233649018930873E-11, 5.26113255735759866E-10, -1.34919775398344888E-09, -1.43632948779513575E-09, 
                    1.03470453927485852E-08, -5.42427480028729824E-09, -4.88475793745928662E-08, 8.40035104689596572E-08, 
                    1.39637754550835535E-07, -5.18773373887414488E-07, -8.48708758607259257E-08, 2.14996026993966525E-06, 
                    -1.55084435011860260E-06, -6.45673042846961899E-06, 1.03735918404559978E-05, 1.33417614992135041E-05, 
                    -4.17514164854039790E-05, -1.15540910383371724E-05, 1.26204335016617082E-04, -4.55568269666841997E-05, 
                    -3.03102046072661175E-04, 2.81763925038067066E-04, 5.81075975053286392E-04, -9.42461407722737748E-04, 
                    -8.44862666553777494E-04, 2.40069778189097322E-03, 7.16982182106401912E-04, -5.07131450921834842E-03, 
                    5.62571574840353152E-04, 9.21478503219718033E-03, -4.13130665603108922E-03, -1.47018820653986824E-02, 
                    1.12904972786859650E-02, 2.09046452556552430E-02, -2.31141340205493172E-02, -2.68914938808945160E-02, 
                    4.00549811051159471E-02, 3.19898775315378067E-02, -6.17662087084131575E-02, -3.66051034028742964E-02, 
                    8.72043982620397534E-02, 4.30958954330476415E-02, -1.14731170710744373E-01, -5.65864586307273792E-02, 
                    1.41414734073382675E-01, 8.56381215561510534E-02, -1.59912565158244369E-01, -1.41795685973059610E-01, 
                    1.49985119618717022E-01, 2.32125963835353111E-01, -6.22665060478243201E-02, -3.21675637808997883E-01, 
                    -1.82867667708335901E-01, 2.13050571355578505E-01, 4.93356078517100782E-01, 4.96591175311718092E-01, 
                    3.30775781411014658E-01, 1.60071993564110698E-01, 5.78899436128592557E-02, 1.56372493475721575E-02, 
                    3.08308811925375173E-03, 4.21170266472711634E-04, 3.57625199426402332E-05, 1.42577664167413176E-06
                };
            }
            default -> throw new IllegalArgumentException("DB" + N + " coefficients not yet implemented");
        };
    }

    /**
     * Normalize coefficients to satisfy Daubechies constraints.
     */
    private static double[] normalizeCoefficients(double[] raw) {
        double tolerance = 1e-10;
        
        // Check if already normalized
        double energy = 0;
        double sum = 0;
        for (double c : raw) {
            energy += c * c;
            sum += c;
        }
        
        // If already normalized (PyWavelets coefficients), return as-is
        if (Math.abs(energy - 1.0) < tolerance && Math.abs(sum - Math.sqrt(2)) < tolerance) {
            return raw.clone();
        }
        
        // Otherwise normalize
        double energyScale = 1.0 / Math.sqrt(energy);
        
        double[] normalized = new double[raw.length];
        for (int i = 0; i < raw.length; i++) {
            normalized[i] = raw[i] * energyScale;
        }
        
        // Then, scale to ensure sum = √2
        sum = 0;
        for (double c : normalized) {
            sum += c;
        }
        double sumScale = Math.sqrt(2) / sum;
        
        for (int i = 0; i < normalized.length; i++) {
            normalized[i] *= sumScale;
        }
        
        return normalized;
    }

    /**
     * Verify that generated coefficients satisfy Daubechies properties.
     * 
     * @param coeffs the coefficients to verify
     * @param N the expected order
     * @return true if coefficients are valid
     * @throws IllegalStateException if a property verification fails with details
     */
    public static boolean verifyDaubechiesProperties(double[] coeffs, int N) {
        if (coeffs.length != 2 * N) {
            throw new IllegalArgumentException(
                String.format("Invalid coefficient length: expected %d, got %d", 2 * N, coeffs.length)
            );
        }
        
        double tolerance = 1e-12;
        
        // 1. Sum of coefficients = √2
        double sum = 0;
        for (double c : coeffs) {
            sum += c;
        }
        double expectedSum = Math.sqrt(2);
        if (Math.abs(sum - expectedSum) > tolerance) {
            throw new IllegalStateException(
                String.format("Sum test failed: %.15f vs %.15f (difference: %.2e, tolerance: %.2e)", 
                    sum, expectedSum, Math.abs(sum - expectedSum), tolerance)
            );
        }
        
        // 2. Sum of squares = 1
        double sumSquares = 0;
        for (double c : coeffs) {
            sumSquares += c * c;
        }
        if (Math.abs(sumSquares - 1.0) > tolerance) {
            throw new IllegalStateException(
                String.format("Sum of squares test failed: %.15f vs 1.0 (difference: %.2e, tolerance: %.2e)", 
                    sumSquares, Math.abs(sumSquares - 1.0), tolerance)
            );
        }
        
        // 3. Orthogonality: sum(h[n] * h[n+2k]) = 0 for k != 0
        for (int k = 1; k < N; k++) {
            double dot = 0;
            for (int n = 0; n < coeffs.length - 2*k; n++) {
                dot += coeffs[n] * coeffs[n + 2*k];
            }
            if (Math.abs(dot) > tolerance) {
                throw new IllegalStateException(
                    String.format("Orthogonality test failed for k=%d: %.15f (tolerance: %.2e)", 
                        k, dot, tolerance)
                );
            }
        }
        
        return true;
    }
}