生物技术

开发平台：

C/C++

blast_seg.c：源码内容

69534.762324, 69543.824745, 69552.887281, 69561.949933, 69571.012701, 69580.075585, 69589.138585, 69598.201700,
69607.264932, 69616.328279, 69625.391742, 69634.455321, 69643.519016, 69652.582827, 69661.646753, 69670.710795,
69679.774953, 69688.839227, 69697.903616, 69706.968121, 69716.032742, 69725.097478, 69734.162330, 69743.227298,
69752.292381, 69761.357580, 69770.422895, 69779.488325, 69788.553871, 69797.619532, 69806.685309, 69815.751202,
69824.817210, 69833.883333, 69842.949572, 69852.015927, 69861.082397, 69870.148982, 69879.215683, 69888.282499,
69897.349431, 69906.416478, 69915.483641, 69924.550919, 69933.618312, 69942.685821, 69951.753445, 69960.821185,
69969.889039, 69978.957009, 69988.025095, 69997.093295, 70006.161611, 70015.230042, 70024.298589, 70033.367250,
70042.436027, 70051.504919, 70060.573926, 70069.643049, 70078.712286, 70087.781639, 70096.851107, 70105.920690,
70114.990388, 70124.060201, 70133.130129, 70142.200172, 70151.270331, 70160.340604, 70169.410992, 70178.481496,
70187.552114, 70196.622848, 70205.693696, 70214.764659, 70223.835738, 70232.906931, 70241.978239, 70251.049662,
70260.121200, 70269.192853, 70278.264621, 70287.336503, 70296.408501, 70305.480613, 70314.552840, 70323.625182,
70332.697639, 70341.770210, 70350.842896, 70359.915697, 70368.988613, 70378.061643, 70387.134788, 70396.208048,
70405.281423, 70414.354912, 70423.428516, 70432.502234, 70441.576067, 70450.650015, 70459.724078, 70468.798254,
70477.872546, 70486.946952, 70496.021473, 70505.096108, 70514.170858, 70523.245722, 70532.320701, 70541.395794,
70550.471002, 70559.546324, 70568.621760, 70577.697311, 70586.772977, 70595.848757, 70604.924651, 70614.000660,
70623.076783, 70632.153020, 70641.229372, 70650.305838, 70659.382418, 70668.459113, 70677.535922, 70686.612845,
70695.689883, 70704.767034, 70713.844300, 70722.921681, 70731.999175, 70741.076784, 70750.154507, 70759.232344,
70768.310295, 70777.388360, 70786.466540, 70795.544833, 70804.623241, 70813.701763, 70822.780399, 70831.859149,
70840.938013, 70850.016991, 70859.096083, 70868.175289, 70877.254609, 70886.334043, 70895.413591, 70904.493253,
70913.573029, 70922.652919, 70931.732923, 70940.813041, 70949.893272, 70958.973618, 70968.054077, 70977.134651,
70986.215338, 70995.296139, 71004.377054, 71013.458082, 71022.539225, 71031.620481, 71040.701851, 71049.783335,
71058.864932, 71067.946644, 71077.028469, 71086.110407, 71095.192460, 71104.274626, 71113.356905, 71122.439299,
71131.521806, 71140.604426, 71149.687160, 71158.770008, 71167.852970, 71176.936045, 71186.019233, 71195.102536,
71204.185951, 71213.269480, 71222.353123, 71231.436879, 71240.520749, 71249.604732, 71258.688829, 71267.773039,
71276.857363, 71285.941799, 71295.026350, 71304.111014, 71313.195791, 71322.280681, 71331.365685, 71340.450802,
71349.536033, 71358.621377, 71367.706834, 71376.792404, 71385.878088, 71394.963885, 71404.049796, 71413.135819,
71422.221956, 71431.308206, 71440.394569, 71449.481045, 71458.567635, 71467.654338, 71476.741154, 71485.828083,
71494.915125, 71504.002280, 71513.089548, 71522.176930, 71531.264424, 71540.352032, 71549.439753, 71558.527586,
71567.615533, 71576.703593, 71585.791766, 71594.880051, 71603.968450, 71613.056962, 71622.145586, 71631.234324,
71640.323174, 71649.412138, 71658.501214, 71667.590403, 71676.679705, 71685.769120, 71694.858648, 71703.948288,
71713.038042, 71722.127908, 71731.217887, 71740.307979, 71749.398184, 71758.488501, 71767.578931, 71776.669474,
71785.760129, 71794.850897, 71803.941778, 71813.032772, 71822.123878, 71831.215097, 71840.306429, 71849.397873,
71858.489430, 71867.581099, 71876.672881, 71885.764776, 71894.856783, 71903.948903, 71913.041135, 71922.133480,
71931.225937, 71940.318507, 71949.411190, 71958.503984, 71967.596892, 71976.689911, 71985.783043, 71994.876288,
72003.969645, 72013.063115, 72022.156696, 72031.250391, 72040.344197, 72049.438116, 72058.532147, 72067.626291,
72076.720547, 72085.814915, 72094.909395, 72104.003988, 72113.098693, 72122.193510, 72131.288440, 72140.383482,
72149.478636, 72158.573902, 72167.669280, 72176.764771, 72185.860373, 72194.956088, 72204.051915, 72213.147854,
72222.243906, 72231.340069, 72240.436344, 72249.532732, 72258.629231, 72267.725843, 72276.822567, 72285.919402,
72295.016350, 72304.113410, 72313.210581, 72322.307865, 72331.405261, 72340.502768, 72349.600388, 72358.698119,
72367.795962, 72376.893918, 72385.991985, 72395.090164, 72404.188455, 72413.286857, 72422.385372, 72431.483998,
72440.582737, 72449.681587, 72458.780548, 72467.879622, 72476.978807, 72486.078104, 72495.177513, 72504.277034,
72513.376666, 72522.476410, 72531.576265, 72540.676233, 72549.776312, 72558.876502, 72567.976805, 72577.077219,
72586.177744, 72595.278381, 72604.379130, 72613.479990, 72622.580962, 72631.682045, 72640.783240, 72649.884547,
72658.985965, 72668.087494, 72677.189135, 72686.290888, 72695.392751, 72704.494727, 72713.596814, 72722.699012,
72731.801321, 72740.903742, 72750.006275, 72759.108919, 72768.211674, 72777.314540, 72786.417518, 72795.520607,
72804.623808, 72813.727120, 72822.830543, 72831.934077, 72841.037723, 72850.141480, 72859.245348, 72868.349327,
72877.453418, 72886.557619, 72895.661932, 72904.766357, 72913.870892, 72922.975538, 72932.080296, 72941.185165,
72950.290145, 72959.395236, 72968.500438, 72977.605751, 72986.711175, 72995.816710, 73004.922356, 73014.028114,
73023.133982, 73032.239962, 73041.346052, 73050.452253, 73059.558566, 73068.664989, 73077.771523, 73086.878168,
73095.984924, 73105.091791, 73114.198769, 73123.305858, 73132.413057, 73141.520368, 73150.627789, 73159.735321,
73168.842964, 73177.950718, 73187.058583, 73196.166558, 73205.274644, 73214.382841, 73223.491149, 73232.599567,
73241.708096, 73250.816736, 73259.925486, 73269.034348, 73278.143320, 73287.252402, 73296.361595, 73305.470899,
73314.580314, 73323.689839, 73332.799474, 73341.909221, 73351.019077, 73360.129045, 73369.239123, 73378.349311,
73387.459610, 73396.570020, 73405.680540, 73414.791170, 73423.901911, 73433.012763, 73442.123725, 73451.234797,
73460.345980, 73469.457273, 73478.568677, 73487.680191, 73496.791815, 73505.903550, 73515.015395, 73524.127351,
73533.239416, 73542.351592, 73551.463879, 73560.576276, 73569.688783, 73578.801400, 73587.914127, 73597.026965,
73606.139913, 73615.252972, 73624.366140, 73633.479419, 73642.592808, 73651.706307, 73660.819916, 73669.933635,
73679.047465, 73688.161404, 73697.275454, 73706.389614, 73715.503884, 73724.618264, 73733.732754, 73742.847354,
73751.962064, 73761.076884, 73770.191814, 73779.306854, 73788.422005, 73797.537265, 73806.652635, 73815.768115,
73824.883705, 73833.999405, 73843.115215, 73852.231135, 73861.347164, 73870.463304, 73879.579553, 73888.695913,
73897.812382, 73906.928961, 73916.045650, 73925.162448, 73934.279357, 73943.396375, 73952.513503, 73961.630741,
73970.748088, 73979.865545, 73988.983112, 73998.100789, 74007.218575, 74016.336472, 74025.454477, 74034.572593,
74043.690818, 74052.809153, 74061.927597, 74071.046151, 74080.164814, 74089.283588, 74098.402470, 74107.521463,
74116.640565, 74125.759776, 74134.879097, 74143.998528, 74153.118068, 74162.237717, 74171.357476, 74180.477344,
74189.597322, 74198.717410, 74207.837607, 74216.957913, 74226.078329, 74235.198854, 74244.319488, 74253.440232,
74262.561085, 74271.682048, 74280.803120, 74289.924301, 74299.045591, 74308.166991, 74317.288500, 74326.410119,
74335.531847, 74344.653684, 74353.775630, 74362.897685, 74372.019850, 74381.142124, 74390.264507, 74399.386999,
74408.509601, 74417.632311, 74426.755131, 74435.878060, 74445.001098, 74454.124245, 74463.247501, 74472.370867,
74481.494341, 74490.617925, 74499.741617, 74508.865419, 74517.989330, 74527.113349, 74536.237478, 74545.361716,
74554.486062, 74563.610518, 74572.735082, 74581.859756, 74590.984538, 74600.109430, 74609.234430, 74618.359539,
74627.484758, 74636.610085, 74645.735520, 74654.861065, 74663.986719, 74673.112481, 74682.238352, 74691.364332,
74700.490421, 74709.616619, 74718.742925, 74727.869340, 74736.995864, 74746.122497, 74755.249238, 74764.376088,
74773.503047, 74782.630114, 74791.757290, 74800.884575, 74810.011969, 74819.139471, 74828.267082, 74837.394801,
74846.522629, 74855.650565, 74864.778611, 74873.906764, 74883.035026, 74892.163397, 74901.291877, 74910.420465,
74919.549161, 74928.677966, 74937.806879, 74946.935901, 74956.065031, 74965.194270, 74974.323617, 74983.453073,
74992.582637, 75001.712309, 75010.842090, 75019.971980, 75029.101977, 75038.232083, 75047.362298, 75056.492620,
75065.623051, 75074.753590, 75083.884238, 75093.014994, 75102.145858, 75111.276831, 75120.407911, 75129.539100,
75138.670397, 75147.801803, 75156.933316, 75166.064938, 75175.196668, 75184.328506, 75193.460453, 75202.592507,
75211.724670, 75220.856940, 75229.989319, 75239.121806, 75248.254401, 75257.387104, 75266.519915, 75275.652835,
75284.785862, 75293.918997, 75303.052240, 75312.185592, 75321.319051, 75330.452618, 75339.586294, 75348.720077,
75357.853968, 75366.987967, 75376.122074, 75385.256289, 75394.390612, 75403.525043, 75412.659582, 75421.794228,
75430.928983, 75440.063845, 75449.198815, 75458.333893, 75467.469078, 75476.604372, 75485.739773, 75494.875282,
75504.010899, 75513.146624, 75522.282456, 75531.418396, 75540.554444, 75549.690599, 75558.826863, 75567.963233,
75577.099712, 75586.236298, 75595.372992, 75604.509793, 75613.646702, 75622.783719, 75631.920843, 75641.058075,
75650.195415, 75659.332862, 75668.470417, 75677.608079, 75686.745848, 75695.883726, 75705.021710, 75714.159802,
75723.298002, 75732.436309, 75741.574724, 75750.713246, 75759.851876, 75768.990613, 75778.129457, 75787.268409,
75796.407468, 75805.546634, 75814.685908, 75823.825290, 75832.964778, 75842.104374, 75851.244077, 75860.383888,
75869.523806, 75878.663831, 75887.803964, 75896.944203, 75906.084550, 75915.225005, 75924.365566, 75933.506235,
75942.647011, 75951.787894, 75960.928884, 75970.069982, 75979.211186, 75988.352498, 75997.493917, 76006.635443,
76015.777076, 76024.918816, 76034.060664, 76043.202618, 76052.344680, 76061.486848, 76070.629124, 76079.771507,
76088.913996, 76098.056593, 76107.199297, 76116.342108, 76125.485025, 76134.628050, 76143.771182, 76152.914420,
76162.057766, 76171.201218, 76180.344777, 76189.488444, 76198.632217, 76207.776097, 76216.920084, 76226.064177,
76235.208378, 76244.352685, 76253.497100, 76262.641621, 76271.786248, 76280.930983, 76290.075824, 76299.220773,
76308.365828, 76317.510989, 76326.656258, 76335.801633, 76344.947114, 76354.092703, 76363.238398, 76372.384200,
76381.530108, 76390.676124, 76399.822245, 76408.968474, 76418.114809, 76427.261251, 76436.407799, 76445.554454,
76454.701215, 76463.848083, 76472.995057, 76482.142138, 76491.289326, 76500.436620, 76509.584021, 76518.731528,
76527.879141, 76537.026861, 76546.174688, 76555.322621, 76564.470660, 76573.618806, 76582.767058, 76591.915417,
76601.063882, 76610.212453, 76619.361131, 76628.509915, 76637.658805, 76646.807802, 76655.956905, 76665.106114,
76674.255430, 76683.404852, 76692.554380, 76701.704015, 76710.853755, 76720.003602, 76729.153556, 76738.303615,
76747.453781, 76756.604052, 76765.754431, 76774.904915, 76784.055505, 76793.206202, 76802.357004, 76811.507913,
76820.658928, 76829.810049, 76838.961276, 76848.112609, 76857.264049, 76866.415594, 76875.567245, 76884.719003,
76893.870866, 76903.022836, 76912.174911, 76921.327092, 76930.479380, 76939.631773, 76948.784273, 76957.936878,
76967.089589, 76976.242406, 76985.395330, 76994.548359, 77003.701493, 77012.854734, 77022.008081, 77031.161533,
77040.315092, 77049.468756, 77058.622526, 77067.776402, 77076.930383, 77086.084471, 77095.238664, 77104.392963,
77113.547368, 77122.701878, 77131.856495, 77141.011216, 77150.166044, 77159.320977, 77168.476017, 77177.631161,
77186.786412, 77195.941768, 77205.097229, 77214.252797, 77223.408470, 77232.564248, 77241.720133, 77250.876122,
77260.032218, 77269.188419, 77278.344725, 77287.501137, 77296.657655, 77305.814278, 77314.971006, 77324.127840,
77333.284780, 77342.441825, 77351.598976, 77360.756232, 77369.913593, 77379.071060, 77388.228632, 77397.386310,
77406.544093, 77415.701982, 77424.859976, 77434.018075, 77443.176279, 77452.334589, 77461.493005, 77470.651525,
77479.810151, 77488.968883, 77498.127719, 77507.286661, 77516.445708, 77525.604860, 77534.764118, 77543.923481,
77553.082949, 77562.242522, 77571.402200, 77580.561984, 77589.721873, 77598.881867, 77608.041966, 77617.202170,
77626.362480, 77635.522894, 77644.683414, 77653.844039, 77663.004769, 77672.165604, 77681.326544, 77690.487589,
77699.648739, 77708.809994, 77717.971354, 77727.132819, 77736.294390, 77745.456065, 77754.617845, 77763.779730,
77772.941720, 77782.103815, 77791.266015, 77800.428320, 77809.590730, 77818.753245, 77827.915864, 77837.078589,
77846.241418, 77855.404353, 77864.567392, 77873.730536, 77882.893784, 77892.057138, 77901.220596, 77910.384160,
77919.547828, 77928.711600, 77937.875478, 77947.039460, 77956.203547, 77965.367739, 77974.532035, 77983.696436,
77992.860942, 78002.025553, 78011.190268, 78020.355088, 78029.520012, 78038.685041, 78047.850175, 78057.015414,
78066.180757, 78075.346204, 78084.511756, 78093.677413, 78102.843174, 78112.009040, 78121.175011, 78130.341086,
78139.507265, 78148.673549, 78157.839938, 78167.006431, 78176.173028, 78185.339730, 78194.506536, 78203.673447,
78212.840462, 78222.007582, 78231.174806, 78240.342135, 78249.509568, 78258.677105, 78267.844746, 78277.012492,
78286.180343, 78295.348297, 78304.516356, 78313.684520, 78322.852787, 78332.021159, 78341.189635, 78350.358216,
78359.526900, 78368.695689, 78377.864583, 78387.033580, 78396.202682, 78405.371887, 78414.541197, 78423.710612,
78432.880130, 78442.049753, 78451.219479, 78460.389310, 78469.559245, 78478.729284, 78487.899427, 78497.069675,
78506.240026, 78515.410481, 78524.581041, 78533.751704, 78542.922472, 78552.093344, 78561.264319, 78570.435399,
78579.606583, 78588.777870, 78597.949262, 78607.120757, 78616.292357, 78625.464060, 78634.635868, 78643.807779,
78652.979794, 78662.151914, 78671.324137, 78680.496464, 78689.668894, 78698.841429, 78708.014068, 78717.186810,
78726.359656, 78735.532606, 78744.705660, 78753.878817, 78763.052079, 78772.225444, 78781.398913, 78790.572485,
78799.746162, 78808.919942, 78818.093826, 78827.267813, 78836.441905, 78845.616100, 78854.790398, 78863.964800,
78873.139306, 78882.313916, 78891.488629, 78900.663446, 78909.838366, 78919.013390, 78928.188518, 78937.363749,
78946.539084, 78955.714522, 78964.890064, 78974.065710, 78983.241458, 78992.417311, 79001.593267, 79010.769326,
79019.945489, 79029.121756, 79038.298125, 79047.474599, 79056.651175, 79065.827856, 79075.004639, 79084.181526,
79093.358517, 79102.535610, 79111.712808, 79120.890108, 79130.067512, 79139.245019, 79148.422630, 79157.600344,
79166.778161, 79175.956081, 79185.134105, 79194.312232, 79203.490462, 79212.668796, 79221.847233, 79231.025773,
79240.204416, 79249.383163, 79258.562012, 79267.740965, 79276.920021, 79286.099181, 79295.278443, 79304.457809,
79313.637277, 79322.816849, 79331.996524, 79341.176302, 79350.356183, 79359.536168, 79368.716255, 79377.896445,
79387.076739, 79396.257135, 79405.437635, 79414.618237, 79423.798943, 79432.979751, 79442.160663, 79451.341678,
79460.522795, 79469.704016, 79478.885339, 79488.066765, 79497.248295, 79506.429927, 79515.611662, 79524.793500,
79533.975441, 79543.157485, 79552.339631, 79561.521881, 79570.704233, 79579.886688, 79589.069246, 79598.251907,
79607.434671, 79616.617537, 79625.800506, 79634.983578, 79644.166753, 79653.350030, 79662.533410, 79671.716893,
79680.900479, 79690.084167, 79699.267958, 79708.451852, 79717.635849, 79726.819948, 79736.004149, 79745.188454,
79754.372861, 79763.557370, 79772.741982, 79781.926697, 79791.111515, 79800.296435, 79809.481457, 79818.666582,
79827.851810, 79837.037140, 79846.222573, 79855.408108, 79864.593746, 79873.779486, 79882.965329, 79892.151274,
79901.337322, 79910.523472, 79919.709725, 79928.896080, 79938.082537, 79947.269097, 79956.455759, 79965.642524,
79974.829391, 79984.016360, 79993.203432, 80002.390606, 80011.577882, 80020.765261, 80029.952742, 80039.140326,
80048.328011, 80057.515799, 80066.703690, 80075.891682, 80085.079777, 80094.267974, 80103.456273, 80112.644675,
80121.833178, 80131.021784, 80140.210492, 80149.399302, 80158.588215, 80167.777229, 80176.966346, 80186.155565,
80195.344886, 80204.534309, 80213.723834, 80222.913462, 80232.103191, 80241.293023, 80250.482956, 80259.672992,
80268.863130, 80278.053369, 80287.243711, 80296.434155, 80305.624700, 80314.815348, 80324.006098, 80333.196950,
80342.387903, 80351.578959, 80360.770116, 80369.961376, 80379.152737, 80388.344201, 80397.535766, 80406.727433,
80415.919202, 80425.111073, 80434.303045, 80443.495120, 80452.687296, 80461.879575, 80471.071955, 80480.264437,
80489.457020, 80498.649706, 80507.842493, 80517.035382, 80526.228373, 80535.421465, 80544.614659, 80553.807955,
80563.001353, 80572.194852, 80581.388453, 80590.582156, 80599.775960, 80608.969866, 80618.163874, 80627.357984,
80636.552195, 80645.746507, 80654.940921, 80664.135437, 80673.330055, 80682.524774, 80691.719594, 80700.914516,
80710.109540, 80719.304665, 80728.499892, 80737.695220, 80746.890650, 80756.086181, 80765.281814, 80774.477548,
80783.673384, 80792.869321, 80802.065360, 80811.261500, 80820.457741, 80829.654084, 80838.850528, 80848.047074,
80857.243721, 80866.440469, 80875.637319, 80884.834270, 80894.031323, 80903.228476, 80912.425732, 80921.623088,
80930.820546, 80940.018105, 80949.215765, 80958.413527, 80967.611390, 80976.809354, 80986.007419, 80995.205586,
81004.403853, 81013.602222, 81022.800693, 81031.999264, 81041.197937, 81050.396710, 81059.595585, 81068.794561,
81077.993638, 81087.192817, 81096.392096, 81105.591477, 81114.790958, 81123.990541, 81133.190225, 81142.390010,
81151.589896, 81160.789883, 81169.989971, 81179.190160, 81188.390450, 81197.590841, 81206.791333, 81215.991926,
81225.192620, 81234.393415, 81243.594311, 81252.795307, 81261.996405, 81271.197604, 81280.398904, 81289.600304,
81298.801805, 81308.003408, 81317.205111, 81326.406915, 81335.608820, 81344.810826, 81354.012932, 81363.215140,
81372.417448, 81381.619857, 81390.822367, 81400.024977, 81409.227689, 81418.430501, 81427.633413, 81436.836427,
81446.039541, 81455.242756, 81464.446072, 81473.649489, 81482.853006, 81492.056624, 81501.260342, 81510.464161,
81519.668081, 81528.872102, 81538.076223, 81547.280444, 81556.484767, 81565.689189, 81574.893713, 81584.098337,
81593.303062, 81602.507887, 81611.712813, 81620.917839, 81630.122966, 81639.328193, 81648.533521, 81657.738949,
81666.944478, 81676.150107, 81685.355837, 81694.561667, 81703.767598, 81712.973629, 81722.179761, 81731.385992,
81740.592325, 81749.798758, 81759.005291, 81768.211924, 81777.418658, 81786.625492, 81795.832427, 81805.039462,
81814.246597, 81823.453833, 81832.661168, 81841.868605, 81851.076141, 81860.283778, 81869.491515, 81878.699352,
81887.907290, 81897.115327, 81906.323465, 81915.531703, 81924.740042, 81933.948480, 81943.157019, 81952.365658,
81961.574397, 81970.783236, 81979.992176, 81989.201215, 81998.410355, 82007.619595, 82016.828935, 82026.038374,
82035.247915, 82044.457555, 82053.667295, 82062.877135, 82072.087075, 82081.297116, 82090.507256, 82099.717496,
82108.927837
};
#define AA20 2
#define LN20 2.9957322735539909
#define CHAR_SET 128
typedef struct SSequence
{
struct SSequence* parent;
char* seq;
Alpha* palpha;
Int4 start;
Int4 length;
Int4 bogus;
Boolean punctuation;
Int4* composition;
Int4* state;
double entropy;
} SSequence;
typedef struct SSeg
{
int begin;
int end;
struct SSeg *next;
} SSeg;
/*---------------------------------------------------------------(SeqNew)---*/
static SSequence* SeqNew(void)
{
SSequence* seq;
seq = (SSequence*) calloc(1, sizeof(SSequence));
if (seq==NULL)
{
/* raise error flag and etc. */
return(seq);
}
seq->parent = (SSequence*) NULL;
seq->seq = (char*) NULL;
seq->palpha = (Alpha*) NULL;
seq->start = seq->length = 0;
seq->bogus = seq->punctuation = FALSE;
seq->composition = seq->state = (Int4*) NULL;
seq->entropy = (double) 0.0;
return(seq);
}
/*------------------------------------------------------------(AlphaFree)---*/
static void AlphaFree (Alpha* palpha)
{
if (!palpha) return;
sfree (palpha->alphaindex);
sfree (palpha->alphaflag);
sfree (palpha->alphachar);
sfree (palpha);
return;
}
/*--------------------------------------------------------------(SeqFree)---*/
static void SeqFree(SSequence* seq)
{
if (seq==NULL) return;
sfree(seq->seq);
AlphaFree(seq->palpha);
sfree(seq->composition);
sfree(seq->state);
sfree(seq);
return;
}
/*--------------------------------------------------------------(SSegFree)---*/
static void SegFree(SSeg* seg)
{
SSeg* nextseg;
while (seg)
{
nextseg = seg->next;
sfree(seg);
seg = nextseg;
}
return;
}
/*--------------------------------------------------------------(hasdash)---*/
static Boolean hasdash(SSequence* win)
{
register char *seq, *seqmax;
seq = win->seq;
seqmax = seq + win->length;
while (seq < seqmax) {
if (*seq++ == '-')
return TRUE;
}
return FALSE;
}
/*------------------------------------------------------------(state_cmp)---*/
static int state_cmp(const void* s1, const void* s2)
{
int *np1, *np2;
np1 = (int *) s1;
np2 = (int *) s2;
return (*np2 - *np1);
}
/*---------------------------------------------------------------(compon)---*/
static void compon(SSequence* win)
{
Int4* comp;
Int4 letter;
char* seq,* seqmax;
Int4* alphaindex;
unsigned char* alphaflag;
Int4 alphasize;
alphasize = win->palpha->alphasize;
alphaindex = win->palpha->alphaindex;
alphaflag = win->palpha->alphaflag;
win->composition = comp =
(Int4*) calloc(alphasize, sizeof(Int4));
seq = win->seq;
seqmax = seq + win->length;
while (seq < seqmax) {
letter = *seq++;
if (!alphaflag[letter])
comp[alphaindex[letter]]++;
else win->bogus++;
}
return;
}
/*--------------------------------------------------------------(stateon)---*/
static void stateon(SSequence* win)
{
Int4 letter, nel, c;
Int4 alphasize;
alphasize = win->palpha->alphasize;
if (win->composition == NULL)
compon(win);
win->state = (Int4*) calloc((alphasize+1), sizeof(win->state[0]));
for (letter = nel = 0; letter < alphasize; ++letter) {
if ((c = win->composition[letter]) == 0)
continue;
win->state[nel++] = c;
}
for (letter = nel; letter < alphasize+1; ++letter)
win->state[letter] = 0;
qsort(win->state, nel, sizeof(win->state[0]), state_cmp);
return;
}
/*--------------------------------------------------------------(openwin)---*/
static SSequence* openwin(SSequence* parent, Int4 start, Int4 length)
{
SSequence* win;
if (start<0 || length<0 || start+length>parent->length)
{
return((SSequence*) NULL);
}
win = (SSequence*) calloc(1, sizeof(SSequence));
/*--- ---[set links, up and down]---*/
win->parent = parent;
win->palpha = parent->palpha;
/*--- ---[install the local copy of the sequence]---*/
win->start = start;
win->length = length;
#if 0 /* Hi Warren! */
win->seq = (char*) calloc(length + 1, sizeof(char));
memcpy(win->seq, (parent->seq)+start, length);
win->seq[length] = '';
#else
win->seq = parent->seq + start;
#endif
win->bogus = 0;
win->punctuation = FALSE;
win->entropy = -2.;
win->state = (Int4*) NULL;
win->composition = (Int4*) NULL;
stateon(win);
return win;
}
/*--------------------------------------------------------------(entropy)---*/
static double entropy(Int4* sv)
{double ent;
Int4 i, total;
total = 0;
for (i=0; sv[i]!=0; i++)
{
total += sv[i];
}
if (total==0) return(0.);
ent = 0.0;
for (i=0; sv[i]!=0; i++)
{
ent += ((double)sv[i])*log(((double)sv[i])/(double)total)/LN2;
}
ent = fabs(ent/(double)total);
return(ent);
}
/*----------------------------------------------------------(decrementsv)---*/
static void decrementsv(Int4* sv, Int4 class)
{
Int4 svi;
while ((svi = *sv++) != 0) {
if (svi == class && *sv < class) {
sv[-1] = svi - 1;
break;
}
}
}
/*----------------------------------------------------------(incrementsv)---*/
static void incrementsv(Int4* sv, Int4 class)
{
for (;;) {
if (*sv++ == class) {
sv[-1]++;
break;
}
}
}
/*------------------------------------------------------------(shiftwin1)---*/
static Int4 shiftwin1(SSequence* win)
{
Int4 j, length;
Int4* comp;
Int4* alphaindex;
unsigned char* alphaflag;
length = win->length;
comp = win->composition;
alphaindex = win->palpha->alphaindex;
alphaflag = win->palpha->alphaflag;
if ((++win->start + length) > win->parent->length) {
--win->start;
return FALSE;
}
if (!alphaflag[j = win->seq[0]])
decrementsv(win->state, comp[alphaindex[j]]--);
else win->bogus--;
j = win->seq[length];
++win->seq;
if (!alphaflag[j])
incrementsv(win->state, comp[alphaindex[j]]++);
else win->bogus++;
if (win->entropy > -2.)
win->entropy = entropy(win->state);
return TRUE;
}
/*-------------------------------------------------------------(closewin)---*/
static void closewin(SSequence* win)
{
if (win==NULL) return;
if (win->state!=NULL) sfree(win->state);
if (win->composition!=NULL) sfree(win->composition);
sfree(win);
return;
}
/*----------------------------------------------------------------(enton)---*/
static void enton(SSequence* win)
{
if (win->state==NULL) {stateon(win);}
win->entropy = entropy(win->state);
return;
}
/*---------------------------------------------------------------(seqent)---*/
static double* seqent(SSequence* seq, Int4 window, Int4 maxbogus)
{
SSequence* win;
double* H;
Int4 i, first, last, downset, upset;
downset = (window+1)/2 - 1;
upset = window - downset;
if (window>seq->length)
{
return((double*) NULL);
}
H = (double*) calloc(seq->length, sizeof(double));
for (i=0; i<seq->length; i++)
{
H[i] = -1.;
}
win = openwin(seq, 0, window);
enton(win);
first = downset;
last = seq->length - upset;
for (i=first; i<=last; i++)
{
if (seq->punctuation && hasdash(win))
{
H[i] = -1.;
shiftwin1(win);
continue;
}
if (win->bogus > maxbogus)
{
H[i] = -1.;
shiftwin1(win);
continue;
}
H[i] = win->entropy;
shiftwin1(win);
}
closewin(win);
return(H);
}
/*------------------------------------------------------------(appendseg)---*/
static void
appendseg(SSeg* segs, SSeg* seg)
{SSeg* temp;
temp = segs;
while (TRUE)
{
if (temp->next==NULL)
{
temp->next = seg;
break;
}
else
{
temp = temp->next;
}
}
return;
}
/*---------------------------------------------------------------(findlo)---*/
static Int4 findlo(Int4 i, Int4 limit, double hicut, double* H)
{
Int4 j;
for (j=i; j>=limit; j--)
{
if (H[j]==-1) break;
if (H[j]>hicut) break;
}
return(j+1);
}
/*---------------------------------------------------------------(findhi)---*/
static Int4 findhi(Int4 i, Int4 limit, double hicut, double* H)
{
Int4 j;
for (j=i; j<=limit; j++)
{
if (H[j]==-1) break;
if (H[j]>hicut) break;
}
return(j-1);
}
/*---------------------------------------------------------------(lnperm)---*/
static double lnperm(Int4* sv, Int4 tot)
{
double ans;
Int4 i;
ans = lnfact[tot];
for (i=0; sv[i]!=0; i++)
{
ans -= lnfact[sv[i]];
}
return(ans);
}
/*----------------------------------------------------------------(lnass)---*/
static double lnass(Int4* sv, Int4 alphasize)
{
double ans;
int svi, svim1;
int class, total;
int i;
ans = lnfact[alphasize];
if (sv[0] == 0)
return ans;
total = alphasize;
class = 1;
svi = *sv;
svim1 = sv[0];
for (i=0;; svim1 = svi) {
if (++i==alphasize) {
ans -= lnfact[class];
break;
}
else if ((svi = *++sv) == svim1) {
class++;
continue;
}
else {
total -= class;
ans -= lnfact[class];
if (svi == 0) {
ans -= lnfact[total];
break;
}
else {
class = 1;
continue;
}
}
}
return ans;
}
/*--------------------------------------------------------------(getprob)---*/
static double getprob(Int4* sv, Int4 total, Alpha* palpha)
{
double ans, ans1, ans2 = 0, totseq;
/* #define LN20 2.9957322735539909 */
/* #define LN4 1.3862943611198906 */
totseq = ((double) total) * palpha->lnalphasize;
/*
ans = lnass(sv, palpha->alphasize) + lnperm(sv, total) - totseq;
*/
ans1 = lnass(sv, palpha->alphasize);
if (ans1 > -100000.0 && sv[0] != INT4_MIN)
{
ans2 = lnperm(sv, total);
}
else
{
#ifdef ERR_POST_EX_DEFINED
ErrPostEx (SEV_ERROR, 0, 0, "Illegal value returned by lnass");
#endif
}
ans = ans1 + ans2 - totseq;
/*
fprintf(stderr, "%lf %lf %lfn", ans, ans1, ans2);
*/
return ans;
}
/*-----------------------------------------------------------------(trim)---*/
static void trim(SSequence* seq, Int4* leftend, Int4* rightend,
SegParameters* sparamsp)
{
SSequence* win;
double prob, minprob;
Int4 shift, len, i;
Int4 lend, rend;
Int4 minlen;
Int4 maxtrim;
/* fprintf(stderr, "%d %dn", *leftend, *rightend); */
lend = 0;
rend = seq->length - 1;
minlen = 1;
maxtrim = sparamsp->maxtrim;
if ((seq->length-maxtrim)>minlen) minlen = seq->length-maxtrim;
minprob = 1.;
/* fprintf(stderr, "n"); -*/
for (len=seq->length; len>minlen; len--)
{
/* fprintf(stderr, "%5d ", len); -*/
win = openwin(seq, 0, len);
i = 0;
shift = TRUE;
while (shift)
{
prob = getprob(win->state, len, sparamsp->palpha);
/* realprob = exp(prob); (for tracing the trim)-*/
/* fprintf(stderr, "%2e ", realprob); -*/
if (prob<minprob)
{
minprob = prob;
lend = i;
rend = len + i - 1;
}
shift = shiftwin1(win);
i++;
}
closewin(win);
/* fprintf(stderr, "n"); -*/
}
/* fprintf(stderr, "%d-%d ", *leftend, *rightend); -*/
*leftend = *leftend + lend;
*rightend = *rightend - (seq->length - rend - 1);
/* fprintf(stderr, "%d-%dn", *leftend, *rightend); -*/
closewin(seq);
return;
}
/*---------------------------------------------------------------(SegSeq)---*/
static void SegSeq(SSequence* seq, SegParameters* sparamsp, SSeg* *segs,
Int4 offset)
{
SSeg* seg,* leftsegs;
SSequence* leftseq;
Int4 window;
double locut, hicut;
Int4 maxbogus;
Int4 downset, upset;
Int4 first, last, lowlim;
Int4 loi, hii, i;
Int4 leftend, rightend, lend, rend;
double* H;
if (sparamsp->window<=0) return;
if (sparamsp->locut<=0.) sparamsp->locut = 0.;
if (sparamsp->hicut<=0.) sparamsp->hicut = 0.;
window = sparamsp->window;
locut = sparamsp->locut;
hicut = sparamsp->hicut;
downset = (window+1)/2 - 1;
upset = window - downset;
maxbogus = sparamsp->maxbogus;
H = seqent(seq, window, maxbogus);
if (H==NULL) return;
first = downset;
last = seq->length - upset;
lowlim = first;
for (i=first; i<=last; i++)
{
if (H[i] <= locut && H[i] != -1.0)
{
loi = findlo(i, lowlim, hicut, H);
hii = findhi(i, last, hicut, H);
leftend = loi - downset;
rightend = hii + upset - 1;
trim(openwin(seq, leftend, rightend-leftend+1), &leftend, &rightend,
sparamsp);
if (i+upset-1<leftend) /* check for trigger window in left trim */
{
lend = loi - downset;
rend = leftend - 1;
leftseq = openwin(seq, lend, rend-lend+1);
leftsegs = (SSeg*) NULL;
SegSeq(leftseq, sparamsp, &leftsegs, offset+lend);
if (leftsegs!=NULL)
{
if (*segs==NULL) *segs = leftsegs;
else appendseg(*segs, leftsegs);
}
closewin(leftseq);
/* trim(openwin(seq, lend, rend-lend+1), &lend, &rend);
seg = (SSeg*) calloc(1, sizeof(SSeg));
seg->begin = lend;
seg->end = rend;
seg->next = (SSeg*) NULL;
if (segs==NULL) segs = seg;
else appendseg(segs, seg); */
}
seg = (SSeg*) calloc(1, sizeof(SSeg));
seg->begin = leftend + offset;
seg->end = rightend + offset;
seg->next = (SSeg*) NULL;
if (*segs==NULL) *segs = seg;
else appendseg(*segs, seg);
i = MIN(hii, rightend+downset);
lowlim = i + 1;
/* i = hii; this ignores the trimmed residues... */
}
}
sfree(H);
return;
}
/*------------------------------------------------------------(mergesegs)---*/
/* merge together overlapping segments,
hilenmin also does something, but we need to ask Scott Federhen what?
*/
static void mergesegs(SSequence* seq, SSeg* segs, Boolean overlaps)
{
SSeg* seg,* nextseg;
Int4 hilenmin; /* hilenmin yet unset */
Int4 len;
hilenmin = 0; /* hilenmin - temporary default */
if (overlaps == FALSE)
return;
if (segs==NULL) return;
if (segs->begin<hilenmin) segs->begin = 0;
seg = segs;
nextseg = seg->next;
while (nextseg!=NULL)
{
if (seg->end>=nextseg->begin && seg->end>=nextseg->end)
{
seg->next = nextseg->next;
sfree(nextseg);
nextseg = seg->next;
continue;
}
if (seg->end>=nextseg->begin) /* overlapping segments */
{
seg->end = nextseg->end;
seg->next = nextseg->next;
sfree(nextseg);
nextseg = seg->next;
continue;
}
len = nextseg->begin - seg->end - 1;
if (len<hilenmin) /* short hient segment */
{
seg->end = nextseg->end;
seg->next = nextseg->next;
sfree(nextseg);
nextseg = seg->next;
continue;
}
seg = nextseg;
nextseg = seg->next;
}
len = seq->length - seg->end - 1;
if (len<hilenmin) seg->end = seq->length - 1;
return;
}
static Int2 SegsToBlastSeqLoc(SSeg* segs, Int4 offset, BlastSeqLoc** seg_locs)
{
SSeqRange* dip;
BlastSeqLoc* last_slp = NULL,* head_slp = NULL;
for ( ; segs; segs = segs->next) {
dip = (SSeqRange*) calloc(1, sizeof(SSeqRange));
dip->left = segs->begin + offset;
dip->right = segs->end + offset;
if (!last_slp) {
last_slp = ListNodeAddPointer(&head_slp, 0, dip);
} else {
last_slp = ListNodeAddPointer(&last_slp, 0, dip);
}
}
*seg_locs = head_slp;
return 0;
}
/*------------------------------------------------------------(AA20alpha)---*/
static Alpha* AA20alpha_std (void)
{
Alpha* palpha;
Int4* alphaindex;
unsigned char* alphaflag;
char* alphachar;
Uint1 c, i;
palpha = (Alpha*) calloc(1, sizeof(Alpha));
palpha->alphabet = AA20;
palpha->alphasize = 20;
palpha->lnalphasize = LN20;
alphaindex = (Int4*) calloc(CHAR_SET , sizeof(Int4));
alphaflag = (unsigned char*) calloc(CHAR_SET , sizeof(char));
alphachar = (char*) calloc(palpha->alphasize , sizeof(char));
for (c=0, i=0; c<128; c++)
{
if (c == 1 || (c >= 3 && c <= 20) || c == 22) {
alphaflag[c] = FALSE;
alphaindex[c] = i;
alphachar[i] = c;
++i;
} else {
alphaflag[c] = TRUE; alphaindex[c] = 20;
}
}
palpha->alphaindex = alphaindex;
palpha->alphaflag = alphaflag;
palpha->alphachar = alphachar;
return (palpha);
}
/*-------------------------------------------------------(SegParametersNewAa)---*/
SegParameters* SegParametersNewAa (void)
{
SegParameters* sparamsp;
sparamsp = (SegParameters*) calloc(1, sizeof(SegParameters));
sparamsp->window = 12;
sparamsp->locut = 2.2;
sparamsp->hicut = 2.5;
sparamsp->period = 1;
sparamsp->hilenmin = 0;
sparamsp->overlaps = FALSE;
sparamsp->maxtrim = 50;
sparamsp->maxbogus = 2;
sparamsp->palpha = AA20alpha_std();
return (sparamsp);
}
/*-------------------------------------------------------(SegParametersCheck)---*/
static void SegParametersCheck (SegParameters* sparamsp)
{
if (!sparamsp) return;
if (sparamsp->window <= 0) sparamsp->window = 12;
if (sparamsp->locut < 0.0) sparamsp->locut = 0.0;
/* if (sparamsp->locut > sparamsp->palpha->lnalphasize)
sparamsp->locut = sparamsp->palpha->lnalphasize; */
if (sparamsp->hicut < 0.0) sparamsp->hicut = 0.0;
/* if (sparamsp->hicut > sparamsp->palpha->lnalphasize)
sparamsp->hicut = sparamsp->palpha->lnalphasize; */
if (sparamsp->locut > sparamsp->hicut)
sparamsp->hicut = sparamsp->locut;
if (sparamsp->maxbogus < 0)
sparamsp->maxbogus = 0;
if (sparamsp->maxbogus > sparamsp->window)
sparamsp->maxbogus = sparamsp->window;
if (sparamsp->period <= 0) sparamsp->period = 1;
if (sparamsp->maxtrim < 0) sparamsp->maxtrim = 0;
return;
}
/*--------------------------------------------------------(SegParametersFree)---*/
void SegParametersFree(SegParameters* sparamsp)
{
if (!sparamsp) return;
AlphaFree(sparamsp->palpha);
sfree(sparamsp);
return;
}
/*------------------------------------------------------------(AlphaCopy)---*/
static Alpha* AlphaCopy (Alpha* palpha)
{
Alpha* pbeta;
Int2 i;
if (!palpha) return((Alpha*) NULL);
pbeta = (Alpha*) calloc(1, sizeof(Alpha));
pbeta->alphabet = palpha->alphabet;
pbeta->alphasize = palpha->alphasize;
pbeta->lnalphasize = palpha->lnalphasize;
pbeta->alphaindex = (Int4*) calloc(CHAR_SET , sizeof(Int4));
pbeta->alphaflag = (unsigned char*) calloc(CHAR_SET , sizeof(char));
pbeta->alphachar = (char*) calloc(palpha->alphasize , sizeof(char));
for (i=0; i<CHAR_SET; i++)
{
pbeta->alphaindex[i] = palpha->alphaindex[i];
pbeta->alphaflag[i] = palpha->alphaflag[i];
}
for (i=0; i<palpha->alphasize; i++)
{
pbeta->alphachar[i] = palpha->alphachar[i];
}
return(pbeta);
}
Int2 SeqBufferSeg (Uint1* sequence, Int4 length, Int4 offset,
SegParameters* sparamsp, BlastSeqLoc** seg_locs)
{
SSequence* seqwin;
SSeg* segs;
Boolean params_allocated = FALSE;
SegParametersCheck (sparamsp);
/* check seg parameters */
if (!sparamsp) {
params_allocated = TRUE;
sparamsp = SegParametersNewAa();
SegParametersCheck (sparamsp);
if (!sparamsp) {
#ifdef ERR_POST_EX_DEFINED
ErrPostEx (SEV_WARNING, 0, 0, "null parameters object");
ErrShow();
#endif
return -1;
}
}
/* make an old-style genwin sequence window object */
seqwin = SeqNew();
seqwin->seq = (char*) sequence;
seqwin->length = length;
seqwin->palpha = AlphaCopy(sparamsp->palpha);
/* seg the sequence */
segs = (SSeg*) NULL;
SegSeq (seqwin, sparamsp, &segs, 0);
/* merge the segment if desired. */
if (sparamsp->overlaps)
mergesegs(seqwin, segs, sparamsp->overlaps);
/* convert segs to seqlocs */
SegsToBlastSeqLoc(segs, offset, seg_locs);
/* clean up & return */
seqwin->seq = NULL;
SeqFree (seqwin);
SegFree (segs);
if(params_allocated)
SegParametersFree(sparamsp);
return 0;
}