/**
 * Determine if atom block can match with the pattern to form a molecule
 * return true if it matches, return false otherwise
 */
static bool try_expand_molecule(t_pack_molecule *molecule,
        const std::multimap<AtomBlockId,t_pack_molecule*>& atom_molecules,
		const AtomBlockId blk_id,
		const t_pack_pattern_block *current_pattern_block) {
	int ipin;
	bool success;
	bool is_optional;
	bool *is_block_optional;
	t_pack_pattern_connections *cur_pack_pattern_connection;
	is_block_optional = molecule->pack_pattern->is_block_optional;
	is_optional = is_block_optional[current_pattern_block->block_id];

    /* If the block in the pattern has already been visited, then there is no need to revisit it */
	if (molecule->atom_block_ids[current_pattern_block->block_id]) {
		if (molecule->atom_block_ids[current_pattern_block->block_id] != blk_id) {
			/* Mismatch between the visited block and the current block implies 
             * that the current netlist structure does not match the expected pattern, 
             * return whether or not this matters */
			return is_optional;
		} else {
			molecule->num_ext_inputs--; /* This block is revisited, implies net is entirely internal to molecule, reduce count */
			return true;
		}
	}

	/* This node has never been visited */
	/* Simplifying assumption: An atom can only map to one molecule */
    auto rng = atom_molecules.equal_range(blk_id);
    bool rng_empty = rng.first == rng.second;
	if(!rng_empty) {
		/* This block is already in a molecule, return whether or not this matters */
		return is_optional;
	}

	if (primitive_type_feasible(blk_id, current_pattern_block->pb_type)) {

		success = true;
		/* If the primitive types match, store it, expand it and explore neighbouring nodes */


        /* store that this node has been visited */
		molecule->atom_block_ids[current_pattern_block->block_id] = blk_id;

		molecule->num_ext_inputs += g_atom_nl.block_input_pins(blk_id).size();
		
		cur_pack_pattern_connection = current_pattern_block->connections;
		while (cur_pack_pattern_connection != NULL && success == true) {
			if (cur_pack_pattern_connection->from_block == current_pattern_block) {
				/* find net corresponding to pattern */
                AtomPortId port_id = g_atom_nl.find_port(blk_id, cur_pack_pattern_connection->from_pin->port->model_port);
                if(!port_id) {
                    //No matching port, we may be at the end
                    success = is_block_optional[cur_pack_pattern_connection->to_block->block_id];
                } else {
                    ipin = cur_pack_pattern_connection->from_pin->pin_number;
                    AtomNetId net_id = g_atom_nl.port_net(port_id, ipin);

                    /* Check if net is valid */
                    if (!net_id || g_atom_nl.net_sinks(net_id).size() != 1) { /* One fanout assumption */
                        success = is_block_optional[cur_pack_pattern_connection->to_block->block_id];
                    } else {
                        auto net_sinks = g_atom_nl.net_sinks(net_id);
                        VTR_ASSERT(net_sinks.size() == 1);

                        auto sink_pin_id = *net_sinks.begin();
                        auto sink_blk_id = g_atom_nl.pin_block(sink_pin_id);

                        success = try_expand_molecule(molecule, atom_molecules, sink_blk_id,
                                    cur_pack_pattern_connection->to_block);
                    }
                }
			} else {
				VTR_ASSERT(cur_pack_pattern_connection->to_block == current_pattern_block);
				/* find net corresponding to pattern */

                auto port_id = g_atom_nl.find_port(blk_id, cur_pack_pattern_connection->to_pin->port->model_port);
                VTR_ASSERT(port_id);
				ipin = cur_pack_pattern_connection->to_pin->pin_number;

                AtomNetId net_id;
				if (cur_pack_pattern_connection->to_pin->port->model_port->is_clock) {
                    VTR_ASSERT(ipin == 1); //TODO: support multi-clock primitives
					net_id = g_atom_nl.port_net(port_id, ipin);
				} else {
					net_id = g_atom_nl.port_net(port_id, ipin);
				}
				/* Check if net is valid */
				if (!net_id || g_atom_nl.net_sinks(net_id).size() != 1) { /* One fanout assumption */
					success = is_block_optional[cur_pack_pattern_connection->from_block->block_id];
				} else {
                    auto driver_blk_id = g_atom_nl.net_driver_block(net_id);
					success = try_expand_molecule(molecule, atom_molecules, driver_blk_id,
                                cur_pack_pattern_connection->from_block);
				}
			}
			cur_pack_pattern_connection = cur_pack_pattern_connection->next;
//.........这里部分代码省略.........

void CGameHelper::GenerateWeaponTargets(const CWeapon* weapon, const CUnit* lastTargetUnit, std::multimap<float, CUnit*>& targets)
{
	GML_RECMUTEX_LOCK(qnum); // GenerateTargets

	const CUnit* attacker = weapon->owner;
	const float radius    = weapon->range;
	const float3& pos     = attacker->pos;
	const float heightMod = weapon->heightMod;
	const float aHeight   = weapon->weaponPos.y;

	// how much damage the weapon deals over 1 second
	const float secDamage = weapon->weaponDef->damages[0] * weapon->salvoSize / weapon->reloadTime * GAME_SPEED;
	const bool paralyzer  = !!weapon->weaponDef->damages.paralyzeDamageTime;

	const std::vector<int>& quads = qf->GetQuads(pos, radius + (aHeight - std::max(0.f, readmap->minheight)) * heightMod);

	const int tempNum = gs->tempNum++;

	typedef std::vector<int>::const_iterator VectorIt;
	typedef std::list<CUnit*>::const_iterator ListIt;
	
	for (VectorIt qi = quads.begin(); qi != quads.end(); ++qi) {
		for (int t = 0; t < teamHandler->ActiveAllyTeams(); ++t) {
			if (teamHandler->Ally(attacker->allyteam, t)) {
				continue;
			}

			const std::list<CUnit*>& allyTeamUnits = qf->GetQuad(*qi).teamUnits[t];

			for (ListIt ui = allyTeamUnits.begin(); ui != allyTeamUnits.end(); ++ui) {
				CUnit* targetUnit = *ui;
				float targetPriority = 1.0f;

				if (luaRules && luaRules->AllowWeaponTarget(attacker->id, targetUnit->id, weapon->weaponNum, weapon->weaponDef->id, &targetPriority)) {
					targets.insert(std::pair<float, CUnit*>(targetPriority, targetUnit));
					continue;
				}


				if (targetUnit->tempNum != tempNum && (targetUnit->category & weapon->onlyTargetCategory)) {
					targetUnit->tempNum = tempNum;

					if (targetUnit->isUnderWater && !weapon->weaponDef->waterweapon) {
						continue;
					}
					if (targetUnit->isDead) {
						continue;
					}

					float3 targPos;
					const unsigned short targetLOSState = targetUnit->losStatus[attacker->allyteam];

					if (targetLOSState & LOS_INLOS) {
						targPos = targetUnit->midPos;
					} else if (targetLOSState & LOS_INRADAR) {
						targPos = targetUnit->midPos + (targetUnit->posErrorVector * radarhandler->radarErrorSize[attacker->allyteam]);
						targetPriority *= 10.0f;
					} else {
						continue;
					}

					const float modRange = radius + (aHeight - targPos.y) * heightMod;

					if ((pos - targPos).SqLength2D() <= modRange * modRange) {
						const float dist2D = (pos - targPos).Length2D();
						const float rangeMul = (dist2D * weapon->weaponDef->proximityPriority + modRange * 0.4f + 100.0f);
						const float damageMul = weapon->weaponDef->damages[targetUnit->armorType] * targetUnit->curArmorMultiple;

						targetPriority *= rangeMul;

						if (targetLOSState & LOS_INLOS) {
							targetPriority *= (secDamage + targetUnit->health);

							if (targetUnit == lastTargetUnit) {
								targetPriority *= weapon->avoidTarget ? 10.0f : 0.4f;
							}

							if (paralyzer && targetUnit->paralyzeDamage > (modInfo.paralyzeOnMaxHealth? targetUnit->maxHealth: targetUnit->health)) {
								targetPriority *= 4.0f;
							}

							if (weapon->hasTargetWeight) {
								targetPriority *= weapon->TargetWeight(targetUnit);
							}
						} else {
							targetPriority *= (secDamage + 10000.0f);
						}

						if (targetLOSState & LOS_PREVLOS) {
							targetPriority /= (damageMul * targetUnit->power * (0.7f + gs->randFloat() * 0.6f));

							if (targetUnit->category & weapon->badTargetCategory) {
								targetPriority *= 100.0f;
							}
							if (targetUnit->crashing) {
								targetPriority *= 1000.0f;
							}
						}

						targets.insert(std::pair<float, CUnit*>(targetPriority, targetUnit));
//.........这里部分代码省略.........

  void AScore::determineHighestScoringPermutations_(const std::vector<std::vector<double> >& peptide_site_scores, std::vector<ProbablePhosphoSites>& sites, const vector<vector<Size> >& permutations, std::multimap<double, Size>& ranking) const
  {
    // For every phospho site of the highest (weighted) scoring phospho site assignment:
    // 1. determine the next best (weighted) score assignment with this site in unphosporylated state.
    // 2. determine the filtering level (peak depths) that maximizes the (unweighted) score difference between these two assignments

    sites.clear();
    // take first set of phospho site assignments
    sites.resize(permutations[0].size());    
    const vector<Size> & best_peptide_sites = permutations[ranking.rbegin()->second]; // sites of the assignment that achieved the highest weighted score

    for (Size i = 0; i < best_peptide_sites.size(); ++i)  // for each phosphorylated site
    {
      multimap<double, Size>::reverse_iterator rev = ranking.rbegin();      
      sites[i].first = best_peptide_sites[i]; // store the site
      sites[i].seq_1 = rev->second; // and permutation
      bool peptide_not_found = true;
      
      // iterate from best scoring peptide to the first peptide that doesn't contain the current phospho site
      do
      {      
        ++rev;
        for (Size j = 0; j < best_peptide_sites.size(); ++j)
        {
          if (j == i)
          {
            if (find(permutations[rev->second].begin(), permutations[rev->second].end(), best_peptide_sites[j]) != permutations[rev->second].end())
            {
              peptide_not_found = true;
              break;
            }
            else
            {
              peptide_not_found = false;
            }
          }
          else
          {
            if (find(permutations[rev->second].begin(), permutations[rev->second].end(), best_peptide_sites[j]) == permutations[rev->second].end())
            {
              peptide_not_found = true;
              break;
            }
            else
            {
              peptide_not_found = false;
            }
          }
        }
      }
      while (peptide_not_found);

      // store permutation of peptide without the phospho site i (seq_2)
      sites[i].seq_2 = rev->second;

      // store phospho site location that is not contained in the best scoring (seq_1) but in seq_2.
      for (Size j = 0; j < permutations[sites[i].seq_2].size(); ++j)
      {
        if (find(permutations[sites[i].seq_1].begin(), permutations[sites[i].seq_1].end(), permutations[sites[i].seq_2][j]) == permutations[sites[i].seq_1].end())
        {
          sites[i].second = permutations[sites[i].seq_2][j];
          break;
        }
      }
    }

    // store peak depth that achieves maximum score difference between best and runner up for every phospho site.
    for (Size i = 0; i < sites.size(); ++i)
    {
      double maximum_score_difference = 0.0;
      sites[i].peak_depth = 1;
      vector<double>::const_iterator first_it = peptide_site_scores[sites[i].seq_1].begin();
      vector<double>::const_iterator second_it = peptide_site_scores[sites[i].seq_2].begin();
      
      for (Size depth = 1; second_it != peptide_site_scores[sites[i].seq_2].end(); ++second_it, ++first_it, ++depth)
      {
        double phospho_at_site_score = *first_it;
        double no_phospho_at_site_score = *second_it;
        double score_difference = phospho_at_site_score - no_phospho_at_site_score;
        
        if (score_difference > maximum_score_difference)
        {
          maximum_score_difference = score_difference;
          sites[i].peak_depth = depth;
        }
      }
    }
  }

bool
MAST::HeatConductionElementBase::
side_external_residual_sensitivity (bool request_jacobian,
                                    RealVectorX& f,
                                    RealMatrixX& jac,
                                    std::multimap<libMesh::boundary_id_type, MAST::BoundaryConditionBase*>& bc) {
    
    
    typedef std::multimap<libMesh::boundary_id_type, MAST::BoundaryConditionBase*> maptype;
    
    // iterate over the boundary ids given in the provided force map
    std::pair<maptype::const_iterator, maptype::const_iterator> it;
    
    const libMesh::BoundaryInfo& binfo = *_system.system().get_mesh().boundary_info;
    
    // for each boundary id, check if any of the sides on the element
    // has the associated boundary
    bool calculate_jac = false;
    
    for (unsigned short int n=0; n<_elem.n_sides(); n++) {
        
        // if no boundary ids have been specified for the side, then
        // move to the next side.
        if (!binfo.n_boundary_ids(&_elem, n))
            continue;
        
        // check to see if any of the specified boundary ids has a boundary
        // condition associated with them
        std::vector<libMesh::boundary_id_type> bc_ids = binfo.boundary_ids(&_elem, n);
        std::vector<libMesh::boundary_id_type>::const_iterator bc_it = bc_ids.begin();
        
        for ( ; bc_it != bc_ids.end(); bc_it++) {
            
            // find the loads on this boundary and evaluate the f and jac
            it = bc.equal_range(*bc_it);
            
            for ( ; it.first != it.second; it.first++) {
                
                // apply all the types of loading
                switch (it.first->second->type()) {
                    case MAST::HEAT_FLUX:
                        calculate_jac = (calculate_jac ||
                                         surface_flux_residual_sensitivity(request_jacobian,
                                                                           f, jac,
                                                                           n,
                                                                           *it.first->second));
                        break;
                        
                    case MAST::CONVECTION_HEAT_FLUX:
                        calculate_jac = (calculate_jac ||
                                         surface_convection_residual_sensitivity(request_jacobian,
                                                                                 f, jac,
                                                                                 n,
                                                                                 *it.first->second));
                        break;
                        
                    case MAST::SURFACE_RADIATION_HEAT_FLUX:
                        calculate_jac = (calculate_jac ||
                                         surface_radiation_residual_sensitivity(request_jacobian,
                                                                                f, jac,
                                                                                n,
                                                                                *it.first->second));
                        break;
                        
                    case MAST::DIRICHLET:
                        // nothing to be done here
                        break;
                        
                    default:
                        // not implemented yet
                        libmesh_error();
                        break;
                }
            }
        }
    }
    return (request_jacobian && calculate_jac);
}

int readFile(const string &filename,
                 std::vector<Interaction> &vecInteractions,
                 std::multimap<size_t, size_t> &mapUserID2Index,
                 std::multimap<size_t, size_t> &mapAdID2Index,
                 std::map<size_t, std::set<size_t> > &mapUserID2AdIDSet){
    //init
    vecInteractions.clear();
    mapUserID2Index.clear();
    mapAdID2Index.clear();

#ifdef DEBUG
    time_t t1, t2;
    t1 = time(NULL);
    std::cout<<"start loading file\n";
    std::cout.flush();
#endif

    //file descriptor
    int fd;
    //file info
    struct stat fileState;
    int status;
    size_t fileSize;

    //the memory-mapped thing itself
    void *mapped;

    //open file for reading
    try{
        fd = open(filename.c_str(), O_RDONLY);
        if(fd < 0){
            throw("failed to open file\n");
        }
    }catch(const char* msg){
        std::cout<<msg<<std::endl;
        std::cout.flush();
        return -1;
    }
    status = fstat(fd, & fileState);
    if(status){
        std::cout<<"failed to get file status\n";

    }
    fileSize = fileState.st_size;

    //memory-map the file
    try{
        //mapped = mmap(0, fileSize, PROT_READ, MAP_PRIVATE| MAP_POPULATE, fd, 0);
        mapped = mmap(0, fileSize, PROT_READ | MAP_POPULATE | MAP_HUGETLB, MAP_SHARED, fd, 0);
        if(mapped == MAP_FAILED){
            throw("failed to map file\n");
        }

    }catch(const char* msg){
        std::cout<<msg<<std::endl;
        std::cout.flush();
        return -1;
    }
    madvise(mapped, fileSize, MADV_SEQUENTIAL|MADV_WILLNEED);




#ifdef DEBUG
    std::cout<<"start conuting line\n";
#endif

    //count line num
    char *ptr = (char*)mapped;
    size_t nChrCnt = 0;
    size_t nLineNum = 0;

    for(; nChrCnt < fileSize; ++nChrCnt){
        if(*ptr == '\n'){
            ++nLineNum;
        }
        ++ptr;
    }

    nChrCnt = 0;
    ptr = static_cast<char*>(mapped);

    vecInteractions.resize(nLineNum);

#ifdef DEBUG
    t2 = time(NULL);
    std::cout << "line num " << nLineNum << " time used " << t2-t1 << "sec \n";
    std::cout.flush();
#endif




    //map userid to its adList
    std::map<size_t, std::set<size_t> >::iterator iter;

    size_t nLineCnt = 0;
    while(1){
        int nEndPos = 0;
        while((nChrCnt < fileSize) && (ptr[nEndPos] != '\n')){
//.........这里部分代码省略.........

void MSFM::StencilS2O2(cv::Mat& image, cv::Mat& u_surface, cv::Mat& state, std::multimap<double,cv::Point>& trial_set, 
        std::map<int,std::multimap<double,cv::Point>::iterator>& mapa_trial, cv::Point& neigh, cv::Point2d &h) {
    
    unsigned char neigh_state = state2M(neigh);
    
    if (neigh_state == P_ALIVE)
        return;
    
    int key = neigh.y + image.rows*neigh.x;
    
    double F = image2M(neigh);
    
    double T_2a = MAX_VAL, T_2b = MAX_VAL, T_1a = MAX_VAL, T_1b = MAX_VAL, gh1, gh2;
    double T_21a = MAX_VAL, T_21b = MAX_VAL, T_11a = MAX_VAL, T_11b = MAX_VAL;
    bool T_1n = false, T_2n = false;
    cv::Point p;
    
    
    gh1 = 1.0/(h.x*h.x + h.y*h.y);
    p = cv::Point(neigh.x - 1, neigh.y - 1);
    if (contains2M(p)) {
        if (state2M(p) == P_ALIVE) {
            T_1a = distance2M(p);
            p = cv::Point(neigh.x - 2, neigh.y - 2);
            if (contains2M(p)) {
                T_11a = distance2M(p);
                if (state2M(p) == P_ALIVE)
                    T_1n = true;
            }
        }
    }
    p = cv::Point(neigh.x + 1, neigh.y + 1);
    if (contains2M(p)) {
        if (state2M(p) == P_ALIVE) {
            T_1b = distance2M(p);
            p = cv::Point(neigh.x + 2, neigh.y + 2);
            if (contains2M(p)) {
                T_11b = distance2M(p);
                if (state2M(p) == P_ALIVE)
                    T_1n = true;
            }
        }
    }
    
    gh2 = 1.0/(h.x*h.x + h.y*h.y);
    p = cv::Point(neigh.x - 1, neigh.y + 1);
    if (contains2M(p)) {
        if (state2M(p) == P_ALIVE) {
            T_2a = distance2M(p);
            p = cv::Point(neigh.x - 2, neigh.y + 2);
            if (contains2M(p)) {
                T_21a = distance2M(p);
                if (state2M(p) == P_ALIVE)
                    T_2n = true;
            }
        }
    }
    p = cv::Point(neigh.x + 1, neigh.y - 1);
    if (contains2M(p)) {
        if (state2M(p) == P_ALIVE) {
            T_2b = distance2M(p);
            p = cv::Point(neigh.x + 2, neigh.y - 2);
            if (contains2M(p)) {
                T_21b = distance2M(p);
                if (state2M(p) == P_ALIVE)
                    T_2n = true;
            }
        }
    }
    
    
    int order = (T_1n && T_2n && ((T_1a >= T_11a && T_1a <= T_1b) || (T_1b >= T_11b && T_1b <= T_1a)) && 
        ((T_2a >= T_21a && T_2a <= T_2b) || (T_2b >= T_21b && T_2b <= T_2a))) ? 2 : 1;
    
    double T_1, T_2;
    switch (order) {
        case 1:
            T_1 = std::min(T_1a, T_1b);
            T_2 = std::min(T_2a, T_2b);
            break;
        case 2:
            if (T_1a >= T_11a && T_1b >= T_11b)
                T_1 = std::min((4.0*T_1a - T_11a)/3.0, (4.0*T_1b - T_11b)/3.0);
            else if (T_1a >= T_11a)
                T_1 = (4.0*T_1a - T_11a)/3.0;
            else
                T_1 = (4.0*T_1b - T_11b)/3.0;
            if (T_2a >= T_21a && T_2b >= T_21b)
                T_2 = std::min((4.0*T_2a - T_21a)/3.0, (4.0*T_2b - T_21b)/3.0);
            else if (T_2a >= T_21a)
                T_2 = (4.0*T_2a - T_21a)/3.0;
            else
                T_2 = (4.0*T_2b - T_21b)/3.0;
            gh1 *= 9.0/4.0;
            gh2 *= 9.0/4.0;
            break;
    }
    
    double cos_phi = fabs(h.x*h.x - h.y*h.y)/(h.x*h.x + h.y*h.y);
    double sin_phi = sqrt(1.0 - cos_phi*cos_phi);
//.........这里部分代码省略.........

/*
* Create an X509_DN
*/
X509_DN::X509_DN(const std::multimap<std::string, std::string>& args)
   {
   for(auto i = args.begin(); i != args.end(); ++i)
      add_attribute(OIDS::lookup(i->first), i->second);
   }

/// ProcessNodesReachableFromGlobals - If we inferred anything about nodes
/// reachable from globals, we have to make sure that we incorporate data for
/// all graphs that include those globals due to the nature of the globals
/// graph.
///
void StructureFieldVisitorBase::
ProcessNodesReachableFromGlobals(DSGraph &DSG,
                                 std::multimap<DSNode*,LatticeValue*> &NodeLVs){
  // Start by marking all nodes reachable from globals.
  DSScalarMap &SM = DSG.getScalarMap();
  if (SM.global_begin() == SM.global_end()) return;

  hash_set<const DSNode*> Reachable;
  for (DSScalarMap::global_iterator GI = SM.global_begin(),
         E = SM.global_end(); GI != E; ++GI)
    SM[*GI].getNode()->markReachableNodes(Reachable);
  if (Reachable.empty()) return;
  
  // If any of the nodes with dataflow facts are reachable from the globals
  // graph, we have to do the GG processing step.
  bool MustProcessThroughGlobalsGraph = false;
  for (std::multimap<DSNode*, LatticeValue*>::iterator I = NodeLVs.begin(),
         E = NodeLVs.end(); I != E; ++I)
    if (Reachable.count(I->first)) {
      MustProcessThroughGlobalsGraph = true;
      break;
    }
  
  if (!MustProcessThroughGlobalsGraph) return;
  Reachable.clear();

  // Compute the mapping from DSG to the globals graph.
  DSGraph::NodeMapTy DSGToGGMap;
  DSG.computeGToGGMapping(DSGToGGMap);

  // Most of the times when we find facts about things reachable from globals we
  // we are in the main graph.  This means that we have *all* of the globals
  // graph in this DSG.  To be efficient, we compute the minimum set of globals
  // that can reach any of the NodeLVs facts.
  //
  // I'm not aware of any wonderful way of computing the set of globals that
  // points to the set of nodes in NodeLVs that is not N^2 in either NodeLVs or
  // the number of globals, except to compute the inverse of DSG.  As such, we
  // compute the inverse graph of DSG, which basically has the edges going from
  // pointed to nodes to pointing nodes.  Because we only care about one
  // connectedness properties, we ignore field info.  In addition, we only
  // compute inverse of the portion of the graph reachable from the globals.
  std::set<std::pair<DSNode*,DSNode*> > InverseGraph;

  for (DSScalarMap::global_iterator GI = SM.global_begin(),
         E = SM.global_end(); GI != E; ++GI)
    ComputeInverseGraphFrom(SM[*GI].getNode(), InverseGraph);

  // Okay, now that we have our bastardized inverse graph, compute the set of
  // globals nodes reachable from our lattice nodes.
  for (std::multimap<DSNode*, LatticeValue*>::iterator I = NodeLVs.begin(),
         E = NodeLVs.end(); I != E; ++I)
    ComputeNodesReachableFrom(I->first, InverseGraph, Reachable);
 
  // Now that we know which nodes point to the data flow facts, figure out which
  // globals point to the data flow facts.
  std::set<GlobalValue*> Globals;
  for (hash_set<const DSNode*>::iterator I = Reachable.begin(),
         E = Reachable.end(); I != E; ++I)
    Globals.insert((*I)->globals_begin(), (*I)->globals_end());

  // Finally, loop over all of the DSGraphs for the program, computing
  // information for the graph if not done already, mapping the result into our
  // context.
  for (hash_map<const Function*, DSGraph*>::iterator GI = ECG.DSInfo.begin(),
         E = ECG.DSInfo.end(); GI != E; ++GI) {
    DSGraph &FG = *GI->second;
    // Graphs can contain multiple functions, only process the graph once.
    if (GI->first != FG.retnodes_begin()->first ||
        // Also, do not bother reprocessing DSG.
        &FG == &DSG)
      continue;

    bool GraphUsesGlobal = false;
    for (std::set<GlobalValue*>::iterator I = Globals.begin(),
           E = Globals.end(); I != E; ++I)
      if (FG.getScalarMap().count(*I)) {
        GraphUsesGlobal = true;
        break;
      }

    // If this graph does not contain the global at all, there is no reason to
    // even think about it.
    if (!GraphUsesGlobal) continue;

    // Otherwise, compute the full set of dataflow effects of the function.
    std::multimap<DSNode*, LatticeValue*> &FGF = getCalleeFacts(FG);
    //std::cerr << "Computed: " << FG.getFunctionNames() << "\n";

#if 0
    for (std::multimap<DSNode*, LatticeValue*>::iterator I = FGF.begin(),
           E = FGF.end(); I != E; ++I)
      I->second->dump();
#endif
    // Compute the mapping of nodes in the globals graph to the function's
//.........这里部分代码省略.........

int main(int, char**)
{
    {
        typedef std::pair<const int, double> V;
        V ar[] =
        {
            V(1, 1),
            V(1, 1.5),
            V(1, 2),
            V(2, 1),
            V(2, 1.5),
            V(2, 2),
            V(3, 1),
            V(3, 1.5),
            V(3, 2),
            V(4, 1),
            V(4, 1.5),
            V(4, 2),
            V(5, 1),
            V(5, 1.5),
            V(5, 2),
            V(6, 1),
            V(6, 1.5),
            V(6, 2),
            V(7, 1),
            V(7, 1.5),
            V(7, 2),
            V(8, 1),
            V(8, 1.5),
            V(8, 2)
        };
        std::multimap<int, double> m(ar, ar+sizeof(ar)/sizeof(ar[0]));
        assert(static_cast<std::size_t>(std::distance(m.begin(), m.end())) == m.size());
        assert(static_cast<std::size_t>(std::distance(m.rbegin(), m.rend())) == m.size());
        std::multimap<int, double>::iterator i;
        i = m.begin();
        std::multimap<int, double>::const_iterator k = i;
        assert(i == k);
        for (int j = 1; j <= 8; ++j)
            for (double d = 1; d <= 2; d += .5, ++i)
            {
                assert(i->first == j);
                assert(i->second == d);
                i->second = 2.5;
                assert(i->second == 2.5);
            }
    }
    {
        typedef std::pair<const int, double> V;
        V ar[] =
        {
            V(1, 1),
            V(1, 1.5),
            V(1, 2),
            V(2, 1),
            V(2, 1.5),
            V(2, 2),
            V(3, 1),
            V(3, 1.5),
            V(3, 2),
            V(4, 1),
            V(4, 1.5),
            V(4, 2),
            V(5, 1),
            V(5, 1.5),
            V(5, 2),
            V(6, 1),
            V(6, 1.5),
            V(6, 2),
            V(7, 1),
            V(7, 1.5),
            V(7, 2),
            V(8, 1),
            V(8, 1.5),
            V(8, 2)
        };
        const std::multimap<int, double> m(ar, ar+sizeof(ar)/sizeof(ar[0]));
        assert(static_cast<std::size_t>(std::distance(m.begin(), m.end())) == m.size());
        assert(static_cast<std::size_t>(std::distance(m.cbegin(), m.cend())) == m.size());
        assert(static_cast<std::size_t>(std::distance(m.rbegin(), m.rend())) == m.size());
        assert(static_cast<std::size_t>(std::distance(m.crbegin(), m.crend())) == m.size());
        std::multimap<int, double>::const_iterator i;
        i = m.begin();
        for (int j = 1; j <= 8; ++j)
            for (double d = 1; d <= 2; d += .5, ++i)
            {
                assert(i->first == j);
                assert(i->second == d);
            }
    }
#if TEST_STD_VER >= 11
    {
        typedef std::pair<const int, double> V;
        V ar[] =
        {
            V(1, 1),
            V(1, 1.5),
            V(1, 2),
            V(2, 1),
            V(2, 1.5),
//.........这里部分代码省略.........

//.........这里部分代码省略.........
    }
    if ((ret = AndroidBitmap_lockPixels(env, bitmapgray, &pixelsgray)) < 0) {
    }
    if ((ret = AndroidBitmap_lockPixels(env, bitmapover, &pixelsover)) < 0) {
    }
    if ((ret = AndroidBitmap_lockPixels(env, bitmapunder, &pixelsunder)) < 0) {
    }
    if ((ret = AndroidBitmap_lockPixels(env, bitmaplaplac, &pixelslaplac)) < 0) {
    }

    //Save
    save_pixelscolor=pixelscolor;
    save_pixelsgray=pixelsgray;
    save_pixelsover=pixelsover;
    save_pixelsunder=pixelsunder;
    save_pixelslaplac=pixelslaplac;

	//pixelscolor=save_pixelscolor;
	graydone(infocolor, pixelscolor, infogray, pixelsgray);
	pixelscolor=save_pixelscolor;
	graydone(infocolor, pixelscolor, infoover, pixelsover);
	pixelscolor=save_pixelscolor;
	graydone(infocolor, pixelscolor, infounder, pixelsunder);
    pixelscolor=save_pixelscolor;
    graydone(infocolor, pixelscolor, infolaplac, pixelslaplac);
    pixelscolor=save_pixelscolor;
    AndroidBitmap_unlockPixels(env, bitmapcolor);

    /*http://imagingsolution.blog107.fc2.com/blog-entry-165.html���sigma��0.8���̗p*/
	double sigma=sgm;
	DifferenceOfGaussianFunc *differenceOfGaussianFunc;
	differenceOfGaussianFunc = new DifferenceOfGaussianFunc(sigma,1);
	pixelsgray=save_pixelsgray;
	gaussian(differenceOfGaussianFunc,infogray,pixelsgray);
	/*Over*/
	sigma=sgm+0.2;
	DifferenceOfGaussianFunc *differenceOfGaussianFuncOver;
	differenceOfGaussianFuncOver = new DifferenceOfGaussianFunc(sigma,1);
	pixelsover=save_pixelsover;
	gaussian(differenceOfGaussianFuncOver,infoover,pixelsover);
	/*Under*/
	sigma=sgm-0.2;
	DifferenceOfGaussianFunc *differenceOfGaussianFuncUnder;
	differenceOfGaussianFuncUnder = new DifferenceOfGaussianFunc(sigma,1);
	pixelsunder=save_pixelsunder;
	gaussian(differenceOfGaussianFuncUnder,infounder,pixelsunder);
	/*Laplac*/
	sigma=sgm;
	DifferenceOfGaussianFunc *differenceOfGaussianFuncLaplac;
	differenceOfGaussianFuncLaplac = new DifferenceOfGaussianFunc(sigma,1);
	pixelslaplac=save_pixelslaplac;
	gaussian(differenceOfGaussianFuncLaplac,infolaplac,pixelslaplac);

    lenamap = (char *)malloc(infogray.height * infogray.width);
    memset(lenamap,'x',infogray.height * infogray.width);
    cntr.clear();
    pixelsgray=save_pixelsgray;

    for (int y=1;y<infogray.height-1;y++) {
        for (int x=1;x<infogray.width-1;x++) {
        	if(strncmp(lenamap+(infogray.width*y)+x, "o", 1) != 0){
        		strncpy(lenamap+(infogray.width*y)+x,"o",1);
        		find_maximal(x,y,pixelsgray,pixelsover,pixelsunder,infogray);
        	}
        }
        pixelsgray = (char *) pixelsgray + infogray.stride;
    }

	/* test */
    pixelsgray=save_pixelsgray;
    std::multimap<int,int>::iterator it = cntr.begin();
    for (int y=1;y<infogray.height-1;y++) {
    	uint8_t * grayline = (uint8_t *) pixelsgray;
    	/*for(int x=1;x<infogray.width-1;x++){
    		grayline[x]=255;
    	}*/
    	//g(differenceOfGaussianFuncUnder,infounder,pixelsunder,(*it).first,(*it).second);
	    	if(y==(*it).first){
		    	while(y==(*it).first){
		    		if(0==g(differenceOfGaussianFuncUnder,infounder,pixelsunder,(*it).first,(*it).second)){
		    			grayline[(*it).second]=255;
		    		}else{
		    			grayline[(*it).second]=0;
		    		}
		    		it++;
		    	}
			}else{
				++it;
			}


		pixelsgray = (char *) pixelsgray + infogray.stride;
    }
    make_orientation(infolaplac,pixelslaplac,cntr);

    AndroidBitmap_unlockPixels(env, bitmapgray);
    AndroidBitmap_unlockPixels(env, bitmapover);
    AndroidBitmap_unlockPixels(env, bitmapunder);
    //free(lenamap);������
}

/// visitGraph - Visit the functions in the specified graph, updating the
/// specified lattice values for all of their uses.
///
void StructureFieldVisitorBase::
visitGraph(DSGraph &DSG, std::multimap<DSNode*, LatticeValue*> &NodeLVs) {
  assert(!NodeLVs.empty() && "No lattice values to compute!");

  // To visit a graph, first step, we visit the instruction making up each
  // function in the graph, but ignore calls when processing them.  We handle
  // call nodes explicitly by looking at call nodes in the graph if needed.  We
  // handle instructions before calls to avoid interprocedural analysis if we
  // can drive lattice values to bottom early.
  //
  SFVInstVisitor IV(DSG, Callbacks, NodeLVs);

  for (DSGraph::retnodes_iterator FI = DSG.retnodes_begin(),
         E = DSG.retnodes_end(); FI != E; ++FI)
    for (Function::iterator BB = FI->first->begin(), E = FI->first->end();
         BB != E; ++BB)
      for (BasicBlock::iterator I = BB->begin(), E = BB->end(); I != E; ++I)
        if (IV.visit(*I) && NodeLVs.empty())
          return;  // Nothing left to analyze.

  // Keep track of which actual direct callees are handled.
  std::set<Function*> CalleesHandled;

  // Once we have visited all of the instructions in the function bodies, if
  // there are lattice values that have not been driven to bottom, see if any of
  // the nodes involved are passed into function calls.  If so, we potentially
  // have to recursively traverse the call graph.
  for (DSGraph::fc_iterator CS = DSG.fc_begin(), E = DSG.fc_end();
       CS != E; ++CS) {
    // Figure out the mapping from a node in the caller (potentially several)
    // nodes in the callee.
    DSGraph::NodeMapTy CallNodeMap;

    Instruction *TheCall = CS->getCallSite().getInstruction();

    // If this is an indirect function call, assume nothing gets passed through
    // it. FIXME: THIS IS BROKEN!  Just get the ECG for the fn ptr if it's not
    // direct.
    if (CS->isIndirectCall())
      continue;

    // If this is an external function call, it cannot be involved with this
    // node, because otherwise the node would be marked incomplete!
    if (CS->getCalleeFunc()->isExternal())
      continue;

    // If we can handle this function call, remove it from the set of direct
    // calls found by the visitor.
    CalleesHandled.insert(CS->getCalleeFunc());

    std::vector<DSNodeHandle> Args;

    DSGraph *CG = &ECG.getDSGraph(*CS->getCalleeFunc());
    CG->getFunctionArgumentsForCall(CS->getCalleeFunc(), Args);

    if (!CS->getRetVal().isNull())
      DSGraph::computeNodeMapping(Args[0], CS->getRetVal(), CallNodeMap);
    for (unsigned i = 0, e = CS->getNumPtrArgs(); i != e; ++i) {
      if (i == Args.size()-1) break;
      DSGraph::computeNodeMapping(Args[i+1], CS->getPtrArg(i), CallNodeMap);
    }
    Args.clear();

    // The mapping we just computed maps from nodes in the callee to nodes in
    // the caller, so we can't query it efficiently.  Instead of going through
    // the trouble of inverting the map to do this (linear time with the size of
    // the mapping), we just do a linear search to see if any affected nodes are
    // passed into this call.
    bool CallCanModifyDataFlow = false;
    for (DSGraph::NodeMapTy::iterator MI = CallNodeMap.begin(),
           E = CallNodeMap.end(); MI != E; ++MI)
      if (NodeLVs.count(MI->second.getNode()))
        // Okay, the node is passed in, check to see if the call might do
        // something interesting to it (i.e. if analyzing the call can produce
        // anything other than "top").
        if ((CallCanModifyDataFlow = NodeCanPossiblyBeInteresting(MI->first,
                                                                  Callbacks)))
          break;

    // If this function call cannot impact the analysis (either because the
    // nodes we are tracking are not passed into the call, or the DSGraph for
    // the callee tells us that analysis of the callee can't provide interesting
    // information), ignore it.
    if (!CallCanModifyDataFlow)
      continue;

    // Okay, either compute analysis results for the callee function, or reuse
    // results previously computed.
    std::multimap<DSNode*, LatticeValue*> &CalleeFacts = getCalleeFacts(*CG);

    // Merge all of the facts for the callee into the facts for the caller.  If
    // this reduces anything in the caller to 'bottom', remove them.
    for (DSGraph::NodeMapTy::iterator MI = CallNodeMap.begin(),
           E = CallNodeMap.end(); MI != E; ++MI) {
      // If we have Lattice facts in the caller for this node in the callee,
      // merge any information from the callee into the caller.

//.........这里部分代码省略.........

bool Polish::Write(const std::string& filename, const std::multimap<int, Airspace>& airspaces) {
	if (airspaces.empty()) {
		AirspaceConverter::LogMessage("Polish output: no airspace, nothing to write", false);
		return false;
	}

	// Check if has the right extension
	if (!boost::iequals(boost::filesystem::path(filename).extension().string(), ".mp")) {
		AirspaceConverter::LogMessage("ERROR: Expected MP extension but found: " + boost::filesystem::path(filename).extension().string(), true);
		return false;
	}

	if (file.is_open()) file.close();
	file.open(filename, std::ios::out | std::ios::trunc | std::ios::binary);
	if (!file.is_open() || file.bad()) {
		AirspaceConverter::LogMessage("ERROR: Unable to open output file: " + filename, true);
		return false;
	}
	AirspaceConverter::LogMessage("Writing output file: " + filename, false);

	WriteHeader(filename);

	// Go trough all airspaces
	for (const std::pair<const int,Airspace>& pair : airspaces)
	{
		// Get the airspace
		const Airspace& a = pair.second;

		// Just a couple if assertions
		assert(a.GetNumberOfPoints() > 3);
		assert(a.GetFirstPoint()==a.GetLastPoint());

		// Determine if it's a POLYGON or a POLYLINE
		if (a.GetType() == Airspace::PROHIBITED || a.GetType() == Airspace::CTR || a.GetType() == Airspace::DANGER) {
			file << "[POLYGON]\n"
				//<< "Type="<< types[a.GetType()] <<"\n"; //TODO...
				<< "Type=0x18" <<"\n";
		} else {
			file << "[POLYLINE]\n"
				<< "Type=0x07\n"; //TODO....
		}

		// Add the label
		file << "Label="<<MakeLabel(a)<<"\n";

		file << "Levels=3\n";

		// Insert all the points
		file << "Data0=";
		double lat,lon;
		for (unsigned int i=0; i<a.GetNumberOfPoints()-1; i++) {
			a.GetPointAt(i).GetLatLon(lat,lon);
			file << "(" << lat << "," << lon << "),";
		}
		a.GetLastPoint().GetLatLon(lat,lon);
		file << "(" << lat << "," << lon << ")\n";

		//file<< "EndLevel=4\n";

		// Close the element
		file << "[END]\n\n";
	}
	file.close();
	return true;
}

/**
 * Return public keys or hashes from scriptPubKey, for 'standard' transaction types.
 */
bool Solver(const CScript& scriptPubKey, txnouttype& typeRet, std::vector<std::vector<unsigned char> >& vSolutionsRet)
{
    // Templates
    static std::multimap<txnouttype, CScript> mTemplates;
    if (mTemplates.empty())
    {
        // Standard tx, sender provides pubkey, receiver adds signature
        mTemplates.insert(std::make_pair(TX_PUBKEY, CScript() << OP_PUBKEY << OP_CHECKSIG));

        // Bitcoin address tx, sender provides hash of pubkey, receiver provides signature and pubkey
        mTemplates.insert(std::make_pair(TX_PUBKEYHASH, CScript() << OP_DUP << OP_HASH160 << OP_PUBKEYHASH << OP_EQUALVERIFY << OP_CHECKSIG));

        // Sender provides N pubkeys, receivers provides M signatures
        mTemplates.insert(std::make_pair(TX_MULTISIG, CScript() << OP_SMALLINTEGER << OP_PUBKEYS << OP_SMALLINTEGER << OP_CHECKMULTISIG));
    }

    vSolutionsRet.clear();

    // Shortcut for pay-to-script-hash, which are more constrained than the other types:
    // it is always OP_HASH160 20 [20 byte hash] OP_EQUAL
    if (scriptPubKey.IsPayToScriptHash())
    {
        typeRet = TX_SCRIPTHASH;
        std::vector<unsigned char> hashBytes(scriptPubKey.begin()+2, scriptPubKey.begin()+22);
        vSolutionsRet.push_back(hashBytes);
        return true;
    }

    int witnessversion;
    std::vector<unsigned char> witnessprogram;
    if (scriptPubKey.IsWitnessProgram(witnessversion, witnessprogram)) {
        if (witnessversion == 0 && witnessprogram.size() == 20) {
            typeRet = TX_WITNESS_V0_KEYHASH;
            vSolutionsRet.push_back(witnessprogram);
            return true;
        }
        if (witnessversion == 0 && witnessprogram.size() == 32) {
            typeRet = TX_WITNESS_V0_SCRIPTHASH;
            vSolutionsRet.push_back(witnessprogram);
            return true;
        }
        return false;
    }

    // Provably prunable, data-carrying output
    //
    // So long as script passes the IsUnspendable() test and all but the first
    // byte passes the IsPushOnly() test we don't care what exactly is in the
    // script.
    if (scriptPubKey.size() >= 1 && scriptPubKey[0] == OP_RETURN && scriptPubKey.IsPushOnly(scriptPubKey.begin()+1)) {
        typeRet = TX_NULL_DATA;
        return true;
    }

    // Scan templates
    const CScript& script1 = scriptPubKey;
    for (const std::pair<txnouttype, CScript>& tplate : mTemplates)
    {
        const CScript& script2 = tplate.second;
        vSolutionsRet.clear();

        opcodetype opcode1, opcode2;
        std::vector<unsigned char> vch1, vch2;

        // Compare
        CScript::const_iterator pc1 = script1.begin();
        CScript::const_iterator pc2 = script2.begin();
        while (true)
        {
            if (pc1 == script1.end() && pc2 == script2.end())
            {
                // Found a match
                typeRet = tplate.first;
                if (typeRet == TX_MULTISIG)
                {
                    // Additional checks for TX_MULTISIG:
                    unsigned char m = vSolutionsRet.front()[0];
                    unsigned char n = vSolutionsRet.back()[0];
                    if (m < 1 || n < 1 || m > n || vSolutionsRet.size()-2 != n)
                        return false;
                }
                return true;
            }
            if (!script1.GetOp(pc1, opcode1, vch1))
                break;
            if (!script2.GetOp(pc2, opcode2, vch2))
                break;

            // Template matching opcodes:
            if (opcode2 == OP_PUBKEYS)
            {
                while (vch1.size() >= 33 && vch1.size() <= 65)
                {
                    vSolutionsRet.push_back(vch1);
                    if (!script1.GetOp(pc1, opcode1, vch1))
                        break;
                }
//.........这里部分代码省略.........