Calculate logo height using the method of:
    
    Schneider TD, Stephens RM. "Sequence logos: a new way to display consensus 
    sequences." Nucleic Acids Res. 1990 Oct 25;18(20):6097-100.
355078
def freqs_to_heights( matrix ):
    
    # Columns are sequence positions, rows are symbol counts/frequencies
    f = matrix.values.transpose()
    n, m = f.shape
    # Ensure normalized
    f = f / sum( f, axis=0 )
    # Shannon entropy (the where replaces 0 with 1 so that '0 log 0 == 0')
    H = - sum( f * log2( where( f, f, 1 ) ), axis=0 )
    # Height
    return transpose( f * ( log2( n ) - H ) )
78217
def height= h
      @head_data.pos = 8
      while bytes = @head_data.read(8)
        length, type = bytes.unpack 'Na*'
        if type == 'IHDR'
          @head_data.pos += 4
          @head_data << [h].pack('N')
          @head_data.pos -= 8
          data = @head_data.read length
          @head_data << [Zlib.crc32(data, Zlib.crc32(type))].pack('N')
          @head_data.rewind
          break
        end
      end
      @head_data.rewind
      h
    end
360321
def height(self, container):
        """"""
        height = float(self.get_style('font_size', container))
        script_level = self.script_level(container)
        if script_level > -1:
            height *= self.position_size * (5 / 6)**script_level
        return height
1335574
public static String toFatCatCore(
			AFPChain afpChain,
			Atom[] ca1,
			Atom[] ca2,
			boolean printLegend, boolean longHeader, boolean showHTML, boolean showAlignmentBlock){

		//TODO The sequence numbers are inaccurate if a !afpChain.isSequential()

		String name1 = afpChain.getName1();
		String name2 = afpChain.getName2();
		int ca1Length = afpChain.getCa1Length();
		int ca2Length = afpChain.getCa2Length();

		int blockNum = afpChain.getBlockNum();
		int totalLenIni = afpChain.getTotalLenIni();
		double totalRmsdIni = afpChain.getTotalRmsdIni();
		int optLength = afpChain.getOptLength();
		double totalRmsdOpt = afpChain.getTotalRmsdOpt();
		double chainRmsd = afpChain.getChainRmsd();
		double alignScore = afpChain.getAlignScore();
		int alnLength = afpChain.getAlnLength();
		int gapLen = afpChain.getGapLen();
		List<AFP> afpSet = afpChain.getAfpSet();

		double similarity = afpChain.getSimilarity();
		double identity = afpChain.getIdentity();

		if (similarity <0  || identity < 0){
			afpChain.calcSimilarity();
			similarity = afpChain.getSimilarity();
			identity = afpChain.getIdentity();
		}



		String algorithmName = afpChain.getAlgorithmName();
		//String version = afpChain.getVersion();

		double probability = afpChain.getProbability();


		int afpNum = afpSet.size();

		int[] blockGap = afpChain.getBlockGap();


		double[] blockScore = afpChain.getBlockScore();
		double[] blockRmsd = afpChain.getBlockRmsd();
		int[] blockSize = afpChain.getBlockSize();


		int alnbeg1 = afpChain.getAlnbeg1();
		int alnbeg2 = afpChain.getAlnbeg2();

		char[] alnseq1 = afpChain.getAlnseq1();
		char[] alnseq2 = afpChain.getAlnseq2();
		char[] alnsymb = afpChain.getAlnsymb();

		// == end of extractation of data values from afpChain
		////////////////////////////////

		StringBuffer txt = new StringBuffer();

		if ( longHeader) {
			txt.append(String.format("Align %s.pdb %d with %s.pdb %d", name1, ca1Length, name2, ca2Length));
		}
		else {
			txt.append(String.format("Align %s.pdb Length1: %d with %s.pdb Length2: %d", name1, ca1Length, name2, ca2Length));
		}
		txt.append(newline);
		if ( afpChain.isShortAlign()){
			txt.append("Short match");
			return txt.toString();
		}
		//txt.append(String.format("raw-score: %.2f norm.-score: %.2f ", alignScore, normAlignScore));

		if ( longHeader ) {
			txt.append(String.format( "Twists %d ini-len %d ini-rmsd %.2f opt-equ %d opt-rmsd %.2f chain-rmsd %.2f Score %.2f align-len %d gaps %d (%.2f%%)",
					blockNum - 1, totalLenIni, totalRmsdIni, optLength, totalRmsdOpt, chainRmsd, alignScore,
					alnLength, gapLen, (100.0 * gapLen/alnLength)) );
			txt.append(newline);

		}  else {

			if ( ! longHeader)
				printScore(txt,algorithmName,probability,longHeader);

			printScoresInLines(afpChain, blockNum, optLength, totalRmsdOpt, alignScore, alnLength, gapLen, identity, similarity,txt);
		}


		//txt.append(String.format("P-value %.2e Afp-num %d Identity %.2f%% Similarity %.2f%% norm.-score: %.2f"+newline, probability, afpNum, identity * 100, similarity * 100, normAlignScore));

		if ( longHeader) {
			printScore(txt,algorithmName,probability,longHeader);

			txt.append(String.format("Afp-num %d Identity %.2f%% Similarity %.2f%%", afpNum, identity * 100, similarity * 100));
			txt.append(newline);
		}

		int i;
		double gap;

		if ( longHeader ){
			int fragLen = 8 ; // FatCatParameters.DEFAULT_FRAGLEN;
			for(i = 0; i < blockNum; i ++)  {
				gap = blockGap[i] /( (double)blockGap[i] + fragLen * blockSize[i]);
				txt.append(String.format( "Block %2d afp %2d score %5.2f rmsd %5.2f gap %d (%.2f%%)",
						i, blockSize[i], blockScore[i], blockRmsd[i], blockGap[i], gap));
				txt.append(newline);
			}
		}

		int     linelen = 70;
		String a;
		String b;
		String c;


		int     t = 0;
		int     ap = alnbeg1;
		int     bp = alnbeg2;
		int     k, len;

		//System.out.println(alnseq1.length + " " + alnseq1.toString());

		while((alnLength - t) > 0)      {
			if(alnLength - t > linelen)     len = linelen;
			else    len = alnLength - t;

			if ( ap >= ca1.length)
				break;
			if ( bp >= ca2.length)
				break;

			String pdb1 = ca1[ap].getGroup().getResidueNumber().toString();
			String pdb2 = ca2[bp].getGroup().getResidueNumber().toString();


			//System.err.println("t,len:"+t+":"+len);
			String lseq1 = new String(alnseq1).substring(t,t+len);
			String lseq2 = new String(alnseq2).substring(t,t+len);
			String lsymb = new String(alnsymb).substring(t,t+len);

			//System.err.println("B:" + b);


			// check conservation and color accordingly, if requested by user.
			if ( showHTML ) {
				a = "";
				b = "";
				c = "";

				//	<span class=\"m\">|</span> ... Structurally equivalent and identical residues
				//  <span class=\"sm\">:</span> ... Structurally equivalent and similar residues
				//  <span class=\"qg\">.</span> ... Structurally equivalent, but not similar residues.

				for (int pos = 0 ; pos < lseq1.length() ; pos ++){
					char c1 = lseq1.charAt(pos);
					char c2 = lseq2.charAt(pos);
					char cl = lsymb.charAt(pos);
					int block = -1 ;
					if ( cl != ' ') {
						try {
							block = Integer.parseInt(String.valueOf(cl));
						} catch (Exception e){
							//
						}
					}
					if ( cl != ' ' ){

						if ( showAlignmentBlock && block > -1 ) {
							a += "<span class=\"alignmentBlock1"+block+"\">" + c1 + "</span>";
							b += "<span class=\"alignmentBlock2"+block+"\">" + c2 + "</span>";
							c += "<span class=\"m\">" + cl + "</span>";
						} else {
							a += getPrefix(c1, c2, 0, block, false).toString() + c1 + "</span>";
							b += getPrefix(c1, c2, 1, block, false).toString() + c2 + "</span>";
							c += "<span class=\"m\">" + cl + "</span>";
						}

					} else if ( c1 != '-' && c2 != '-') {

						a += "<span class=\"sm\">" + c1 + "</span>";
						b += "<span class=\"sm\">" + c2 + "</span>";
						c += "<span class=\"sm\">" + cl + "</span>";


					} else {

						a += "<span class=\"qg\">" + c1 + "</span>";
						b += "<span class=\"qg\">" + c2 + "</span>";
						c += "<span class=\"qg\">" + cl + "</span>";

					}

					if(c1 != '-') ap ++;
					if(c2 != '-') bp ++;
				}


			} else {

				a = lseq1;
				b = lseq2;
				c = lsymb;
			}

			txt.append(newline);
			if ( longHeader )
				txt.append(String.format("%14s", " "));
			else
				txt.append(String.format("%14s", " "));

			if (  longHeader ) {
				for(k = 10; k <= len; k += 10)
					txt.append("    .    :");
				if(k <= len + 5) txt.append("    .");

			} else {

				for(k = 10; k <= len; k += 10)
					txt.append("----+----|");
				if(k <= len + 5) txt.append("----+");


			}



			txt.append(newline);
			txt.append(String.format("Chain 1:%5s %s"+newline +"%14s%s"+newline+"Chain 2:%5s %s",
					pdb1, a, " ", c, pdb2, b));

			txt.append(newline);

			if ( ! showHTML){
				for(k = 0; k < len; k ++)       {
					if(a.charAt(k) != '-') ap ++;
					if(b.charAt(k) != '-') bp ++;
				}
			}
			t += len;



		}
		txt.append(newline);
		if ( printLegend ){
			if ( algorithmName.equalsIgnoreCase(CeMain.algorithmName) ||
					algorithmName.equalsIgnoreCase(SmithWaterman3Daligner.algorithmName)){
				txt.append("Note: positions are from PDB; | means alignment of identical amino acids, : of similar amino acids ");

			} else {
				txt.append("Note: positions are from PDB; the numbers between alignments are block index");
			}
			txt.append(newline);

		}
		return txt.toString();

	}

	public static void printScoresInLines(AFPChain afpChain, int blockNum, int optLength, double totalRmsdOpt, double alignScore,
			int alnLength, int gapLen, double identity, double similarity, StringBuffer txt)
	{
		if ( blockNum - 1 > 0) {
			txt.append(String.format( "Twists %d ", blockNum -1 ));
			txt.append(newline);
		}

		txt.append(String.format("Equ: %d ", optLength));
		txt.append(newline);
		txt.append(String.format("RMSD: %.2f ", totalRmsdOpt));
		txt.append(newline);
		txt.append(String.format("Score: %.2f ", alignScore));
		txt.append(newline);
		txt.append(String.format("Align-len: %d ", alnLength));
		txt.append(newline);
		txt.append(String.format("Gaps: %d (%.2f%%)",
				gapLen, (100.0 * gapLen/alnLength)) );
		txt.append(newline);
		if ( afpChain.getTMScore() >= 0) {
			txt.append(String.format("TM-score: %.2f",afpChain.getTMScore()));
			txt.append(newline);
		}
		txt.append(newline);
		txt.append(String.format("Identity: %.2f%% ", identity * 100 ));
		txt.append(newline);
		txt.append(String.format("Similarity: %.2f%%", similarity * 100));
		txt.append(newline);
	}

	private static void printScore(StringBuffer txt,
			String algorithmName,
			double probability,
			boolean longHeader)
	{
		if ( algorithmName.equalsIgnoreCase(CeMain.algorithmName) || algorithmName.equalsIgnoreCase(CeSideChainMain.algorithmName) ){
			txt.append(String.format("Z-score %.2f ", probability));
			if ( ! longHeader)
				txt.append(newline);
		} else if ( algorithmName.equalsIgnoreCase(SmithWaterman3Daligner.algorithmName)) {

		} else {
			if ( longHeader ){
				txt.append(String.format("P-value %.2e ",probability));
			}  else {
				txt.append(String.format("P-value: %.2e ",probability));
				txt.append(newline);
			}
		}



	}

	/**
	 * Prints the afpChain as a nicely formatted alignment, including alignment
	 * statistics, the aligned sequences themselves, and information about the
	 * superposition.
	 * @param afpChain
	 * @param ca1
	 * @param ca2
	 * @return a String representation as it is used on the RCSB PDB web site for display.
	 */
	public static String toWebSiteDisplay(AFPChain afpChain, Atom[] ca1, Atom[] ca2 ){
		boolean showAlignmentBlock   = true;
		return toWebSiteDisplay(afpChain, ca1, ca2, showAlignmentBlock);
	}

	/**
	 * Prints the afpChain as a nicely formatted alignment, including alignment
	 * statistics, the aligned sequences themselves, and information about the
	 * superposition.
	 * @param afpChain
	 * @param ca1
	 * @param ca2
	 *
	 * @return a String representation as it is used on the RCSB PDB web site for display.
	 */
	public static String toWebSiteDisplay(AFPChain afpChain, Atom[] ca1, Atom[] ca2, boolean showAlignmentBlock){

		boolean printLegend = true;
		boolean longHeader  = true;
		boolean showHTML = true;

		if ( afpChain.getAlgorithmName().equalsIgnoreCase(FatCatFlexible.algorithmName)) {


			String msg =  toFatCatCore(afpChain,ca1,ca2,printLegend,longHeader,showHTML, showAlignmentBlock);

			return msg;
		}

		boolean showSeq = true;

		AFPAlignmentDisplay.getAlign(afpChain, ca1, ca2, showSeq);


		//      String msg= toFatCatCore(afpChain,ca1,ca2, printLegend,longHeader);
		//


		String msg = toPrettyAlignment(afpChain, ca1, ca2, showHTML, showAlignmentBlock);


		msg = msg + newline +
		"     <span class=\"m\">|</span> ... Structurally equivalent and identical residues " + newline +
		"     <span class=\"sm\">:</span> ... Structurally equivalent and similar residues " + newline +
		"     <span class=\"qg\">.</span> ... Structurally equivalent, but not similar residues. " + newline;

		msg += newline;
		msg += "     To calculate the coordinates of chain 2 aligned on chain 1 apply the following transformation: ";
		msg += newline;
		msg += newline;
		msg += toRotMat(afpChain);
		return msg;


	}

	public static String toPrettyAlignment(AFPChain afpChain, Atom[] ca1, Atom[] ca2, boolean showHTML, boolean showAlignmentBlock) {
		String name1 = afpChain.getName1();
		String name2 = afpChain.getName2();
		int ca1Length = afpChain.getCa1Length();
		int ca2Length = afpChain.getCa2Length();

		int blockNum = afpChain.getBlockNum();


		int optLength = afpChain.getOptLength();
		double totalRmsdOpt = afpChain.getTotalRmsdOpt();

		double alignScore = afpChain.getAlignScore();
		int alnLength = afpChain.getAlnLength();
		int gapLen = afpChain.getGapLen();


		double similarity = afpChain.getSimilarity();
		double identity = afpChain.getIdentity();

		if (similarity <0  || identity < 0){
			afpChain.calcSimilarity();
			similarity = afpChain.getSimilarity();
			identity = afpChain.getIdentity();
		}


		String algorithmName = afpChain.getAlgorithmName();
		//String version = afpChain.getVersion();

		double probability = afpChain.getProbability();


		// == end of extractation of data values from afpChain

		StringBuffer txt = new StringBuffer();

		txt.append(String.format("Align %s.pdb Length1: %d with %s.pdb Length2: %d", name1, ca1Length, name2, ca2Length));

		txt.append(newline);

		if ( afpChain.isShortAlign()){
			txt.append("Short match");
			return txt.toString();
		}

		printScore(txt, algorithmName, probability, false);
		printScoresInLines(afpChain, blockNum, optLength, totalRmsdOpt, alignScore, alnLength, gapLen,identity, similarity, txt);
		txt.append(newline);

		int[] optLen = afpChain.getOptLen();
		int[][][] optAln = afpChain.getOptAln();


		int i, j,p1, p2;

		int k;
		int p1b = 0;
		int p2b = 0;

		int     len = 0;
		StringWriter alnseq1 = new StringWriter();
		StringWriter alnseq2 = new StringWriter();
		StringWriter alnsymb = new StringWriter();
		StringWriter header1  = new StringWriter();
		StringWriter footer1  = new StringWriter();
		StringWriter header2  = new StringWriter();
		StringWriter footer2  = new StringWriter();
		StringWriter block    = new StringWriter();

		int aligPos = -1;
		for(i = 0; i < blockNum; i ++)  {

			for(j = 0; j < optLen[i]; j ++) {

				p1 = optAln[i][0][j];
				p2 = optAln[i][1][j];
				aligPos++;

				//               System.out.println(p1 + " " + p2 + " " +  footer2.toString());

				if ( len == 0){
					//the first position of sequence in alignment
					formatStartingText(p1,p2,header1,header2,footer1,footer2,ca1,ca2);
				} else {
					// check for gapped region
					int lmax = (p1 - p1b - 1)>(p2 - p2b - 1)?(p1 - p1b - 1):(p2 - p2b - 1);
					for(k = 0; k < lmax; k ++)      {


						formatGappedRegion(ca1, ca2, txt, p1, p2, k, p1b, p2b, alnseq1, alnseq2, alnsymb, header1, footer1, header2,
								footer2, block,len, showHTML);
						len++;
						doLenCheck(len,txt,header1,header2,alnseq1,alnsymb,alnseq2,footer1, footer2,block, showHTML)  ;
					}
				}

				// ALIGNED REGION
				//           System.out.println(len + " >" + header1.toString() +"< ");
				//           System.out.println(len + " >" + header2.toString() +"< ");
				//           System.out.println(len + " >" + alnseq1.toString() +"< ");
				//           System.out.println(len + " >" + alnsymb.toString() +"< ");
				//           System.out.println(len + " >" + alnseq2.toString() +"< ");
				//           System.out.println(len + " >" + footer1.toString() +"< ");
				formatAlignedRegion(afpChain, ca1, ca2, p1, p2, alnseq1, alnseq2, alnsymb, header1, footer1, header2, footer2, block,len, aligPos, showHTML, showAlignmentBlock);
				//            System.out.println(len + " >" + header1.toString() +"< ");
				//            System.out.println(len + " >" + header2.toString() +"< ");
				//            System.out.println(len + " >" + alnseq1.toString() +"< ");
				//            System.out.println(len + " >" + alnsymb.toString() +"< ");
				//            System.out.println(len + " >" + alnseq2.toString() +"< ");
				//            System.out.println(len + " >" + footer1.toString() +"< ");

				len++;

				doLenCheck(len,txt,header1,header2,alnseq1,alnsymb,alnseq2,footer1, footer2,block, showHTML)  ;

				p1b = p1;
				p2b = p2;

				//header1.append(newline);
				//header2.append(newline);

			}

		}

		alnLength = len;

		doLenCheck(LINELENGTH,txt,header1,header2,alnseq1,alnsymb,alnseq2,footer1, footer2,block, showHTML);
		return txt.toString();
	}

	/**
	 * Prints the alignment in the simplest form: a list of aligned residues.
	 * Format is one line per residue pair, tab delimited:
	 * <ul><li>1. PDB number. Includes insertion code</li>
	 * <li>1. Chain.</li>
	 * <li>1. Amino Acid. Three letter code.</li>
	 * <li>2. PDB number.</li>
	 * <li>2. Chain.</li>
	 * <li>2. Amino Acid.</li>
	 * </ul>
	 * example:
	 * <code>152	A	ALA	161S	A	VAL</code>
	 * <p>Note that this format loses information about blocks.
	 * @param afpChain
	 * @param ca1
	 * @param ca2
	 * @return a String representation of the aligned pairs.
	 */
	public static String toAlignedPairs(AFPChain afpChain, Atom[] ca1, Atom[] ca2) {
		StringWriter pairs = new StringWriter();

		//Write structure names & PDB codes
		pairs.append("#Struct1:\t");
		pairs.append(afpChain.getName1());
		pairs.append("\n");
		pairs.append("#Struct2:\t");
		pairs.append(afpChain.getName2());
		pairs.append("\n");

		//Write optimally aligned pairs
		pairs.append("#Num1\tChain1\tAA1\tNum2\tChain2\tAA2\n");
		int[][][] optAln = afpChain.getOptAln();
		int[] blockLen = afpChain.getOptLen();
		for( int block=0;block<afpChain.getBlockNum(); block++) {
			for(int i=0;i<blockLen[block];i++) {
				Atom atom1 = ca1[ optAln[block][0][i] ];
				Atom atom2 = ca2[ optAln[block][1][i] ];

				pairs.append(atom1.getGroup().getResidueNumber().toString());
				pairs.append('\t');
				pairs.append(atom1.getGroup().getChain().getName());
				pairs.append('\t');
				pairs.append(atom1.getGroup().getPDBName());
				pairs.append('\t');
				pairs.append(atom2.getGroup().getResidueNumber().toString());
				pairs.append('\t');
				pairs.append(atom2.getGroup().getChain().getName());
				pairs.append('\t');
				pairs.append(atom2.getGroup().getPDBName());
				pairs.append('\n');
			}
		}

		return pairs.toString();
	}

	private static void formatGappedRegion(Atom[] ca1, Atom[] ca2, StringBuffer txt, int p1, int p2, int k, int p1b, int p2b,
			StringWriter alnseq1, StringWriter alnseq2, StringWriter alnsymb, StringWriter header1, StringWriter footer1,
			StringWriter header2, StringWriter footer2, StringWriter block, int len, boolean formatHTML)	{

		// DEAL WITH GAPS
		int tmppos = (p1 - p1b - 1);
		block.append("g");

		int  pos1=p1b+1+k ;
		char oneletter1 = ' ';
		try {
			oneletter1 = getOneLetter(ca1[pos1].getGroup());
		} catch (Exception e){}
		int pos2=p2b+1+k;
		char oneletter2 = ' ';
		try {
			oneletter2 = getOneLetter(ca2[pos2].getGroup());
		} catch (Exception e){}


		if(k >= tmppos) {
			//alnseq1[len] = '-';
			if (  formatHTML){
				alnseq1.append("<span class=\"qg\">-</span>");
				header1.append(" ");
				header2.append(" ");


			} else {
				alnseq1.append('-');
				header1.append(" ");
				header2.append(" ");
			}

		}
		else {
			if ( formatHTML){
				alnseq1.append(getPrefix(oneletter1,oneletter2,0,-1, false));
			}
			alnseq1.append(oneletter1);
			if (formatHTML){
				alnseq1.append("</span>");
			}
			formatPosition(pos1,ca1, len, header1, header2);

		}

		if(k >= (p2 - p2b - 1)) {
			//alnseq2[len] = '-';
			if ( formatHTML){
				alnseq2.append("<span class=\"qg\">-</span>");
				footer1.append(" ");
				footer2.append(" ");
			} else {
				alnseq2.append('-');
				footer1.append(" ");
				footer2.append(" ");
			}

		}
		else  {
			if ( formatHTML){
				alnseq2.append(getPrefix(oneletter1,oneletter2,1, -1, false));
			}
			alnseq2.append(oneletter2);
			if (formatHTML){
				alnseq2.append("</span>");
			}
			formatPosition(pos2, ca2, len, footer1, footer2);

		}
		//alnsymb[len ++] = ' ';
		alnsymb.append(' ');

	}
774988
def find_consensus(bases):
    
    nucs = ['A', 'T', 'G', 'C', 'N']
    total = sum([bases[nuc] for nuc in nucs if nuc in bases])
    # save most common base as consensus (random nuc if there is a tie)
    try:
        top = max([bases[nuc] for nuc in nucs if nuc in bases])
    except:
        bases['consensus'] = ('N', 'n/a')
        bases['consensus frequency'] = 'n/a'
        bases['reference frequency'] = 'n/a'
        return bases
    top = [(nuc, bases[nuc]) for nuc in bases if bases[nuc] == top]
    if top[0][1] == 0:
        bases['consensus'] = ('n/a', 0)
    else:
        bases['consensus'] = random.choice(top)
    if total == 0:
        c_freq = 'n/a'
        ref_freq = 'n/a'
    else:
        c_freq = float(bases['consensus'][1]) / float(total)
        if bases['ref'] not in bases:
            ref_freq = 0
        else:
            ref_freq = float(bases[bases['ref']]) / float(total)
    bases['consensus frequency'] = c_freq
    bases['reference frequency'] = ref_freq
    return bases
689308
def get_block_height(height, api_code=None):
    

    resource = 'block-height/{0}?format=json'.format(height)
    if api_code is not None:
        resource += '&api_code=' + api_code
    response = util.call_api(resource)
    json_response = json.loads(response)
    return [Block(b) for b in json_response['blocks']]
1350255
private static Rectangle calcRotatedSize(int width, int height, int degree) {
		if (degree >= 90) {
			if (degree / 90 % 2 == 1) {
				int temp = height;
				height = width;
				width = temp;
			}
			degree = degree % 90;
		}
		double r = Math.sqrt(height * height + width * width) / 2;
		double len = 2 * Math.sin(Math.toRadians(degree) / 2) * r;
		double angel_alpha = (Math.PI - Math.toRadians(degree)) / 2;
		double angel_dalta_width = Math.atan((double) height / width);
		double angel_dalta_height = Math.atan((double) width / height);
		int len_dalta_width = (int) (len * Math.cos(Math.PI - angel_alpha - angel_dalta_width));
		int len_dalta_height = (int) (len * Math.cos(Math.PI - angel_alpha - angel_dalta_height));
		int des_width = width + len_dalta_width * 2;
		int des_height = height + len_dalta_height * 2;

		return new Rectangle(des_width, des_height);
	}
179872
def find_n50(self):
        
        for sample in self.metadata:
            # Initialise the N50 attribute in case there is no assembly, and the attribute is not created in the loop
            sample[self.analysistype].n50 = '-'
            # Initialise a variable to store a running total of contig lengths
            currentlength = 0
            for contig_length in sample[self.analysistype].contig_lengths:
                # Increment the current length with the length of the current contig
                currentlength += contig_length
                # If the current length is now greater than the total genome / 2, the current contig length is the N50
                if currentlength >= sample[self.analysistype].genome_length * 0.5:
                    # Populate the dictionary, and break the loop
                    sample[self.analysistype].n50 = contig_length
                    break
1272196
def get_conservation(block):
    
    consensus = block['sequences'][0]['seq']
    assert all(c.isupper() for c in consensus), \
            "So-called consensus contains indels!"
    # remove all non-consensus positions -- now alignment is easy
    cleaned = [[c for c in s['seq'] if not c.islower()]
               for s in block['sequences'][1:]]
    height = float(len(cleaned))
    # validation
    for row in cleaned:
        if len(row) != len(consensus):
            raise ValueError("Aligned sequence length (%s) doesn't match "
                             "consensus (%s)"
                             % (len(row), len(consensus)))
    # transpose & go
    columns = zip(*cleaned)
    return dict((idx + 1, columns[idx].count(cons_char) / height)
                for idx, cons_char in enumerate(consensus))
909165
function(sequence, n, startSymbol, endSymbol, stats) {
    var result = [];
    frequencies = {};
    nrOfNgrams = 0;
    
    if (!_(sequence).isArray()) {
        sequence = tokenizer.tokenize(sequence);
    }

    var count = _.max([0, sequence.length - n + 1]);

    // Check for left padding    
    if(typeof startSymbol !== "undefined" && startSymbol !== null) {
        // Create an array of (n) start symbols
        var blanks = [];
        for(var i = 0 ; i < n ; i++) {
            blanks.push(startSymbol);
        }

        // Create the left padding
        for(var p = n - 1 ; p > 0 ; p--) {
            // Create a tuple of (p) start symbols and (n - p) words
          var ngram = blanks.slice(0, p).concat(sequence.slice(0, n - p));
          result.push(ngram);
          if (stats) {
            countNgrams(ngram);
          }
        }
    }

    // Build the complete ngrams
    for (var i = 0; i < count; i++) {
        var ngram = sequence.slice(i, i + n);
        result.push(ngram);
        if (stats) { 
          countNgrams(ngram);
        }
    }

    // Check for right padding
    if(typeof endSymbol !== "undefined" && endSymbol !== null) {
        // Create an array of (n) end symbols
        var blanks = [];
        for(var i = 0 ; i < n ; i++) {
            blanks.push(endSymbol);
        }

        // create the right padding
        for(var p = n - 1 ; p > 0 ; p--) {
            // Create a tuple of (p) start symbols and (n - p) words
          var ngram = sequence.slice(sequence.length - p, sequence.length).concat(blanks.slice(0, n - p));
          result.push(ngram);
          if (stats) {
            countNgrams(ngram);
          }
        }
    }
    
    if (stats) {
      
      // Count frequencies
      var Nr = {};
      Object.keys(frequencies).forEach(function(key) {
        if (!Nr[frequencies[key]]) {
          Nr[frequencies[key]] = 0;
        }
        Nr[frequencies[key]]++;
      });
      
      // Return the ngrams AND statistics
      return {
        "ngrams": result,
        "frequencies": frequencies,
        "Nr": Nr,
        "numberOfNgrams": nrOfNgrams
      };
      
    }
    else { // Do not break existing API of this module 
      return result;
    }
}
135750
def print_consensus(genomes):
    
    # generate consensus sequences
    cons = {} # cons[genome][sample][contig] = consensus
    for genome, contigs in list(genomes.items()):
        cons[genome] = {}
        for contig, samples in list(contigs.items()):
            for sample, stats in list(samples.items()):
                if sample not in cons[genome]:
                    cons[genome][sample] = {}
                seq = cons[genome][sample][contig] = []
                for pos, ps in enumerate(stats['bp_stats'], 1):
                    ref, consensus = ps['ref'], ps['consensus'][0]
                    if consensus == 'n/a':
                        consensus = ref.lower()
                    seq.append(consensus)
    # print consensus sequences
    for genome, samples in cons.items():
        for sample, contigs in samples.items():
            fn = '%s.%s.consensus.fa' % (genome, sample)
            f = open(fn, 'w')
            for contig, seq in contigs.items():
                print('>%s' % (contig), file = f)
                print(''.join(seq), file = f)
            f.close()
    return cons
177476
def calculate_dimensions(image, long_side, short_side):
    """"""
    if image.width >= image.height:
        return '{0}x{1}'.format(long_side, short_side)
    return '{0}x{1}'.format(short_side, long_side)
1033488
def get_symbol(num_classes, num_layers, image_shape, num_group=32, conv_workspace=256, dtype='float32', **kwargs):
    
    image_shape = [int(l) for l in image_shape.split(',')]
    (nchannel, height, width) = image_shape
    if height <= 32:
        num_stages = 3
        if (num_layers-2) % 9 == 0 and num_layers >= 164:
            per_unit = [(num_layers-2)//9]
            filter_list = [16, 64, 128, 256]
            bottle_neck = True
        elif (num_layers-2) % 6 == 0 and num_layers < 164:
            per_unit = [(num_layers-2)//6]
            filter_list = [16, 16, 32, 64]
            bottle_neck = False
        else:
            raise ValueError("no experiments done on num_layers {}, you can do it yourself".format(num_layers))
        units = per_unit * num_stages
    else:
        if num_layers >= 50:
            filter_list = [64, 256, 512, 1024, 2048]
            bottle_neck = True
        else:
            filter_list = [64, 64, 128, 256, 512]
            bottle_neck = False
        num_stages = 4
        if num_layers == 18:
            units = [2, 2, 2, 2]
        elif num_layers == 34:
            units = [3, 4, 6, 3]
        elif num_layers == 50:
            units = [3, 4, 6, 3]
        elif num_layers == 101:
            units = [3, 4, 23, 3]
        elif num_layers == 152:
            units = [3, 8, 36, 3]
        elif num_layers == 200:
            units = [3, 24, 36, 3]
        elif num_layers == 269:
            units = [3, 30, 48, 8]
        else:
            raise ValueError("no experiments done on num_layers {}, you can do it yourself".format(num_layers))

    return resnext(units      = units,
                  num_stages  = num_stages,
                  filter_list = filter_list,
                  num_classes = num_classes,
                  num_group   = num_group,
                  image_shape = image_shape,
                  bottle_neck = bottle_neck,
                  workspace   = conv_workspace,
                  dtype       = dtype)
1140752
def get_orfs(_orf_length = 100, prediction = @prediction)
      '-' if prediction.type != 'nucleotide'

      seq = Bio::Sequence::NA.new(prediction.raw_sequence)

      result = {}
      key = 0
      (1..6).each do |f|
        s = seq.translate(f)
        f = -1 if f == 4
        f = -2 if f == 5
        f = -3 if f == 6
        s.scan(/(\w{30,})/) do |_orf|
          orf_start = $LAST_MATCH_INFO.offset(0)[0] + 1
          orf_end   = $LAST_MATCH_INFO.offset(0)[1] + 1
          coverage = (((orf_end - orf_start) / s.length.to_f) * 100).ceil
          # reduce the orf_end and the translated length by 2% to increase the
          #   width between ORFs on the plot
          chopping = s.length * 0.02
          orf_end = (orf_end.to_f - chopping).floor
          translated_length = (s.length - chopping).ceil
          key += 1
          result[key] = { frame: f, orf_start: orf_start, orf_end: orf_end,
                          coverage: coverage,
                          translated_length: translated_length }
        end
      end
      result
    end
599392
def _read_species_and_inital_geometry(self):
        
        header_pattern = r"Standard Nuclear Orientation \(Angstroms\)\s+I\s+Atom\s+X\s+Y\s+Z\s+-+"
        table_pattern = r"\s*\d+\s+([a-zA-Z]+)\s*([\d\-\.]+)\s*([\d\-\.]+)\s*([\d\-\.]+)\s*"
        footer_pattern = r"\s*-+"
        temp_geom = read_table_pattern(self.text, header_pattern,
                                       table_pattern, footer_pattern)
        if temp_geom == None or len(temp_geom) == 0:
            self.data["species"] = None
            self.data["initial_geometry"] = None
            self.data["initial_molecule"] = None
            self.data["point_group"] = None
        else:
            temp_point_group = read_pattern(
            self.text, {
                "key":
                r"Molecular Point Group\s+([A-Za-z\d\*]+)"
            },
            terminate_on_match=True).get('key')
            if temp_point_group != None:
                self.data["point_group"] = temp_point_group[0][0]
            else:
                self.data["point_group"] = None
            temp_geom = temp_geom[0]
            species = []
            geometry = np.zeros(shape=(len(temp_geom), 3), dtype=float)
            for ii, entry in enumerate(temp_geom):
                species += [entry[0]]
                for jj in range(3):
                    geometry[ii, jj] = float(entry[jj + 1])
            self.data["species"] = species
            self.data["initial_geometry"] = geometry
            self.data["initial_molecule"] = Molecule(
                species=species,
                coords=geometry,
                charge=self.data.get('charge'),
                spin_multiplicity=self.data.get('multiplicity'))
960637
protected function _calculate_height()
    {
        $style = $this->_frame->get_style();
        $height = $style->height;

        $cellmap = $this->_frame->get_cellmap();
        $cellmap->assign_frame_heights();
        $rows = $cellmap->get_rows();

        // Determine our content height
        $content_height = 0;
        foreach ($rows as $r) {
            $content_height += $r["height"];
        }

        $cb = $this->_frame->get_containing_block();

        if (!($style->overflow === "visible" ||
            ($style->overflow === "hidden" && $height === "auto"))
        ) {
            // Only handle min/max height if the height is independent of the frame's content

            $min_height = $style->min_height;
            $max_height = $style->max_height;

            if (isset($cb["h"])) {
                $min_height = $style->length_in_pt($min_height, $cb["h"]);
                $max_height = $style->length_in_pt($max_height, $cb["h"]);

            } else if (isset($cb["w"])) {
                if (mb_strpos($min_height, "%") !== false) {
                    $min_height = 0;
                } else {
                    $min_height = $style->length_in_pt($min_height, $cb["w"]);
                }
                if (mb_strpos($max_height, "%") !== false) {
                    $max_height = "none";
                } else {
                    $max_height = $style->length_in_pt($max_height, $cb["w"]);
                }
            }

            if ($max_height !== "none" && $min_height > $max_height) {
                // Swap 'em
                list($max_height, $min_height) = array($min_height, $max_height);
            }

            if ($max_height !== "none" && $height > $max_height) {
                $height = $max_height;
            }

            if ($height < $min_height) {
                $height = $min_height;
            }
        } else {
            // Use the content height or the height value, whichever is greater
            if ($height !== "auto") {
                $height = $style->length_in_pt($height, $cb["h"]);

                if ($height <= $content_height) {
                    $height = $content_height;
                } else {
                    $cellmap->set_frame_heights($height, $content_height);
                }
            } else {
                $height = $content_height;
            }
        }

        return $height;
    }
705491
def _my_eps_formatter(logodata, format, ordered_alphabets) :
    
    substitutions = {}
    from_format =[
        "creation_date",    "logo_width",           "logo_height",      
        "lines_per_logo",   "line_width",           "line_height",
        "line_margin_right","line_margin_left",     "line_margin_bottom",
        "line_margin_top",  "title_height",         "xaxis_label_height",
        "creator_text",     "logo_title",           "logo_margin",
        "stroke_width",     "tic_length",           
        "stacks_per_line",  "stack_margin",
        "yaxis_label",      "yaxis_tic_interval",   "yaxis_minor_tic_interval",
        "xaxis_label",      "xaxis_tic_interval",   "number_interval",
        "fineprint",        "shrink_fraction",      "errorbar_fraction",
        "errorbar_width_fraction",
        "errorbar_gray",    "small_fontsize",       "fontsize",
        "title_fontsize",   "number_fontsize",      "text_font",
        "logo_font",        "title_font",          
        "logo_label",       "yaxis_scale",          "end_type",
        "debug",            "show_title",           "show_xaxis",
        "show_xaxis_label", "show_yaxis",           "show_yaxis_label",
        "show_boxes",       "show_errorbars",       "show_fineprint",
        "rotate_numbers",   "show_ends",            "stack_height",
        "stack_width"
        ]
   
    for s in from_format :
        substitutions[s] = getattr(format,s)

    substitutions["shrink"] = str(format.show_boxes).lower()


    # --------- COLORS --------------
    def format_color(color):
        return  " ".join( ("[",str(color.red) , str(color.green), 
            str(color.blue), "]"))  

    substitutions["default_color"] = format_color(format.default_color)

    colors = []  
    if hasattr(format.color_scheme, 'rules'):
        grouplist = format.color_scheme.rules
    else:
        # this line needed for weblogo 3.4
        grouplist = format.color_scheme.groups
    for group in grouplist:
        cf = format_color(group.color)
        for s in group.symbols :
            colors.append( "  ("+s+") " + cf )
    substitutions["color_dict"] = "\n".join(colors)
        
    data = []
    
    # Unit conversion. 'None' for probability units
    conv_factor = None #JDB
    #JDB conv_factor = std_units[format.unit_name]
    
    data.append("StartLine")

    seq_from = format.logo_start- format.first_index
    seq_to = format.logo_end - format.first_index +1

    # seq_index : zero based index into sequence data
    # logo_index : User visible coordinate, first_index based
    # stack_index : zero based index of visible stacks
    for seq_index in range(seq_from, seq_to) :
        logo_index = seq_index + format.first_index 
        stack_index = seq_index - seq_from
        
        if stack_index!=0 and (stack_index % format.stacks_per_line) ==0 :
            data.append("")
            data.append("EndLine")
            data.append("StartLine")
            data.append("")
        
        data.append("(%s) StartStack" % format.annotate[seq_index] )

        if conv_factor: 
            stack_height = logodata.entropy[seq_index] * std_units[format.unit_name]
        else :
            stack_height = 1.0 # Probability

        # The following code modified by JDB to use ordered_alphabets
        # and also to replace the "blank" characters 'b' and 'B'
        # by spaces.
        s_d = dict(zip(logodata.alphabet, logodata.counts[seq_index]))
        s = []
        for aa in ordered_alphabets[seq_index]:
            if aa not in ['B', 'b']:
                s.append((s_d[aa], aa))
            else:
                s.append((s_d[aa], ' '))
#        s = [(s_d[aa], aa) for aa in ordered_alphabets[seq_index]]

        # Sort by frequency. If equal frequency then reverse alphabetic
        # (So sort reverse alphabetic first, then frequencty)
        # TODO: doublecheck this actual works
        #s = list(zip(logodata.counts[seq_index], logodata.alphabet))
        #s.sort(key= lambda x: x[1])
        #s.reverse()
        #s.sort(key= lambda x: x[0])
        #if not format.reverse_stacks: s.reverse()

        C = float(sum(logodata.counts[seq_index])) 
        if C > 0.0 :
            fraction_width = 1.0
            if format.scale_width :
                fraction_width = logodata.weight[seq_index] 
            # print(fraction_width, file=sys.stderr)
            for c in s:
                data.append(" %f %f (%s) ShowSymbol" % (fraction_width, c[0]*stack_height/C, c[1]) )

        # Draw error bar on top of logo. Replaced by DrawErrorbarFirst above.
        if logodata.entropy_interval is not None and conv_factor and C>0.0:

            low, high = logodata.entropy_interval[seq_index]
            center = logodata.entropy[seq_index]
            low *= conv_factor
            high *= conv_factor
            center *=conv_factor
            if high> format.yaxis_scale : high = format.yaxis_scale 

            down = (center - low) 
            up   = (high - center) 
            data.append(" %f %f DrawErrorbar" % (down, up) )
            
        data.append("EndStack")
        data.append("")
               
    data.append("EndLine")
    substitutions["logo_data"] = "\n".join(data)  


    # Create and output logo
    template = corebio.utils.resource_string( __name__, '_weblogo_template.eps', __file__).decode()
    logo = string.Template(template).substitute(substitutions)

    return logo.encode()