C++ GeoExtent::getSRS方法代码示例

osg::Node*
SingleKeyNodeFactory::createTile(TileModel* model, bool setupChildrenIfNecessary)
{
    // compile the model into a node:
    TileNode* tileNode = _modelCompiler->compile( model, _frame );

    // see if this tile might have children.
    bool prepareForChildren =
        setupChildrenIfNecessary &&
        model->_tileKey.getLOD() < *_options.maxLOD();

    osg::Node* result = 0L;

    if ( prepareForChildren )
    {
        //Compute the min range based on the 2D size of the tile
        osg::BoundingSphere bs = tileNode->getBound();
        GeoExtent extent = model->_tileKey.getExtent();
        GeoPoint lowerLeft(extent.getSRS(), extent.xMin(), extent.yMin(), 0.0, ALTMODE_ABSOLUTE);
        GeoPoint upperRight(extent.getSRS(), extent.xMax(), extent.yMax(), 0.0, ALTMODE_ABSOLUTE);
        osg::Vec3d ll, ur;
        lowerLeft.toWorld( ll );
        upperRight.toWorld( ur );
        double radius = (ur - ll).length() / 2.0;
        float minRange = (float)(radius * _options.minTileRangeFactor().value());

        TilePagedLOD* plod = new TilePagedLOD( _engineUID, _liveTiles, _deadTiles );
        plod->setCenter  ( bs.center() );
        plod->addChild   ( tileNode );
        plod->setRange   ( 0, minRange, FLT_MAX );
        plod->setFileName( 1, Stringify() << tileNode->getKey().str() << "." << _engineUID << ".osgearth_engine_mp_tile" );
        plod->setRange   ( 1, 0, minRange );

#if USE_FILELOCATIONCALLBACK
        osgDB::Options* options = Registry::instance()->cloneOrCreateOptions();
        options->setFileLocationCallback( new FileLocationCallback() );
        plod->setDatabaseOptions( options );
#endif
        
        result = plod;

        // this one rejects back-facing tiles:
        if ( _frame.getMapInfo().isGeocentric() && _options.clusterCulling() == true )
        {
            osg::HeightField* hf =
                model->_elevationData.getHeightField();

            result->addCullCallback( HeightFieldUtils::createClusterCullingCallback(
                hf,
                tileNode->getKey().getProfile()->getSRS()->getEllipsoid(),
                *_options.verticalScale() ) );
        }
    }
    else
    {
        result = tileNode;
    }

    return result;
}

osg::BoundingSphered
FeatureModelGraph::getBoundInWorldCoords( const GeoExtent& extent ) const
{
    osg::Vec3d center, corner;

    GeoExtent workingExtent;

    if ( extent.getSRS()->isEquivalentTo( _usableMapExtent.getSRS() ) )
    {
        workingExtent = extent;
    }
    else
    {
        workingExtent = extent.transform( _usableMapExtent.getSRS() ); // safe.
    }

    workingExtent.getCentroid( center.x(), center.y() );
    corner.x() = workingExtent.xMin();
    corner.y() = workingExtent.yMin();

    if ( _session->getMapInfo().isGeocentric() )
    {
        workingExtent.getSRS()->transformToECEF( center, center );
        workingExtent.getSRS()->transformToECEF( corner, corner );
    }

    return osg::BoundingSphered( center, (center-corner).length() );
}

void
UnifiedCubeProfile::getIntersectingTiles(const GeoExtent&      remoteExtent,
                                         unsigned              localLOD,
                                         std::vector<TileKey>& out_intersectingKeys ) const
{
    if ( getSRS()->isHorizEquivalentTo( remoteExtent.getSRS() ) )
    {
        addIntersectingTiles( remoteExtent, localLOD, out_intersectingKeys );
    }
    else
    {
        // the cube profile is non-contiguous. so there may be multiple local extents required
        // to fully intersect the remote extent.

        // first transform the remote extent to lat/long.
        GeoExtent remoteExtent_gcs = remoteExtent.getSRS()->isGeographic()
            ? remoteExtent
            : remoteExtent.transform( remoteExtent.getSRS()->getGeographicSRS() );

        // Chop the input extent into three separate extents: for the equatorial, north polar,
        // and south polar tile regions.
        for( int face=0; face<6; ++face )
        {
            GeoExtent partExtent_gcs = _faceExtent_gcs[face].intersectionSameSRS( remoteExtent_gcs );
            if ( partExtent_gcs.isValid() )
            {
                GeoExtent partExtent = transformGcsExtentOnFace( partExtent_gcs, face );
                addIntersectingTiles( partExtent, localLOD, out_intersectingKeys );
            }
        }        
    }
}

osg::Node*
SerialKeyNodeFactory::createTile(TileModel* model,
                                 bool       tileHasRealData)
{
    // compile the model into a node:
    TileNode* tileNode = _modelCompiler->compile( model );

    // see if this tile might have children.
    bool prepareForChildren =
        (tileHasRealData || (_options.minLOD().isSet() && model->_tileKey.getLOD() < *_options.minLOD())) &&
        model->_tileKey.getLOD() < *_options.maxLOD();

    osg::Node* result = 0L;

    if ( prepareForChildren )
    {
        //Compute the min range based on the 2D size of the tile
        osg::BoundingSphere bs = tileNode->getBound();
        GeoExtent extent = model->_tileKey.getExtent();
        GeoPoint lowerLeft(extent.getSRS(), extent.xMin(), extent.yMin(), 0.0, ALTMODE_ABSOLUTE);
        GeoPoint upperRight(extent.getSRS(), extent.xMax(), extent.yMax(), 0.0, ALTMODE_ABSOLUTE);
        osg::Vec3d ll, ur;
        lowerLeft.toWorld( ll );
        upperRight.toWorld( ur );
        double radius = (ur - ll).length() / 2.0;
        float minRange = (float)(radius * _options.minTileRangeFactor().value());

        osgDB::Options* dbOptions = Registry::instance()->cloneOrCreateOptions();

        TileGroup* plod = new TileGroup(tileNode, _engineUID, _liveTiles.get(), _deadTiles.get(), dbOptions);
        plod->setSubtileRange( minRange );


#if USE_FILELOCATIONCALLBACK
        dbOptions->setFileLocationCallback( new FileLocationCallback() );
#endif
        
        result = plod;
    }
    else
    {
        result = tileNode;
    }

    // this one rejects back-facing tiles:
    if ( _mapInfo.isGeocentric() && _options.clusterCulling() == true )
    {
        osg::HeightField* hf =
            model->_elevationData.getHeightField();

        result->addCullCallback( HeightFieldUtils::createClusterCullingCallback(
            hf,
            tileNode->getKey().getProfile()->getSRS()->getEllipsoid(),
            *_options.verticalScale() ) );
    }

    return result;
}

osg::BoundingSphered
FeatureModelGraph::getBoundInWorldCoords(const GeoExtent& extent,
                                         const MapFrame*  mapf ) const
{
    osg::Vec3d center, corner;
    GeoExtent workingExtent;

    if ( !extent.isValid() )
    {
        return osg::BoundingSphered();
    }

    if ( extent.getSRS()->isEquivalentTo( _usableMapExtent.getSRS() ) )
    {
        workingExtent = extent;
    }
    else
    {
        workingExtent = extent.transform( _usableMapExtent.getSRS() ); // safe.
    }

    workingExtent.getCentroid( center.x(), center.y() );
    
    double centerZ = 0.0;    
    if ( mapf )
    {
        // Use an appropriate resolution for this extents width
        double resolution = workingExtent.width();
        ElevationQuery query( *mapf );
        GeoPoint p( mapf->getProfile()->getSRS(), center, ALTMODE_ABSOLUTE );
        query.getElevation( p, center.z(), resolution );
        centerZ = center.z();
    }    

    corner.x() = workingExtent.xMin();
    corner.y() = workingExtent.yMin();
    corner.z() = 0;

    if ( _session->getMapInfo().isGeocentric() )
    {
        const SpatialReference* ecefSRS = workingExtent.getSRS()->getECEF();
        workingExtent.getSRS()->transform( center, ecefSRS, center );
        workingExtent.getSRS()->transform( corner, ecefSRS, corner );
        //workingExtent.getSRS()->transformToECEF( center, center );
        //workingExtent.getSRS()->transformToECEF( corner, corner );
    }

    if (workingExtent.getSRS()->isGeographic() &&
        ( workingExtent.width() >= 90 || workingExtent.height() >= 90 ) )
    {
        return osg::BoundingSphered( osg::Vec3d(0,0,0), 2*center.length() );
    }

    return osg::BoundingSphered( center, (center-corner).length() );
}

void
Profile::getIntersectingTiles(const GeoExtent& extent, unsigned localLOD, std::vector<TileKey>& out_intersectingKeys) const
{
    GeoExtent ext = extent;

    // reproject into the profile's SRS if necessary:
    if ( !getSRS()->isHorizEquivalentTo( extent.getSRS() ) )
    {
        // localize the extents and clamp them to legal values
        ext = clampAndTransformExtent( extent );
        if ( !ext.isValid() )
            return;
    }

    if ( ext.crossesAntimeridian() )
    {
        GeoExtent first, second;
        if (ext.splitAcrossAntimeridian( first, second ))
        {
            addIntersectingTiles( first, localLOD, out_intersectingKeys );
            addIntersectingTiles( second, localLOD, out_intersectingKeys );
        }
    }
    else
    {
        addIntersectingTiles( ext, localLOD, out_intersectingKeys );
    }
}

osg::BoundingSphere SimplePager::getBounds(const TileKey& key) const
{
    int samples = 6;

    GeoExtent extent = key.getExtent();

    double xSample = extent.width() / (double)samples;
    double ySample = extent.height() / (double)samples;

    osg::BoundingSphere bs;
    for (int c = 0; c < samples+1; c++)
    {
        double x = extent.xMin() + (double)c * xSample;
        for (int r = 0; r < samples+1; r++)
        {
            double y = extent.yMin() + (double)r * ySample;
            osg::Vec3d world;

            GeoPoint samplePoint(extent.getSRS(), x, y, 0, ALTMODE_ABSOLUTE);

            GeoPoint wgs84 = samplePoint.transform(osgEarth::SpatialReference::create("epsg:4326"));
            wgs84.toWorld(world);
            bs.expandBy(world);
        }
    }
    return bs;
}

bool
VerticalDatum::transform(const VerticalDatum* from,
                         const VerticalDatum* to,
                         const GeoExtent&     extent,
                         osg::HeightField*    hf )
{
    if ( from == to )
        return true;

    unsigned cols = hf->getNumColumns();
    unsigned rows = hf->getNumRows();
    
    osg::Vec3d sw(extent.west(), extent.south(), 0.0);
    osg::Vec3d ne(extent.east(), extent.north(), 0.0);
    
    double xstep = abs(extent.east() - extent.west()) / double(cols-1);
    double ystep = abs(extent.north() - extent.south()) / double(rows-1);

    if ( !extent.getSRS()->isGeographic() )
    {
        const SpatialReference* geoSRS = extent.getSRS()->getGeographicSRS();
        extent.getSRS()->transform(sw, geoSRS, sw);
        extent.getSRS()->transform(ne, geoSRS, ne);
        xstep = (ne.x()-sw.x()) / double(cols-1);
        ystep = (ne.y()-sw.y()) / double(rows-1);
    }

    for( unsigned c=0; c<cols; ++c)
    {
        double lon = sw.x() + xstep*double(c);
        for( unsigned r=0; r<rows; ++r)
        {
            double lat = sw.y() + ystep*double(r);
            float& h = hf->getHeight(c, r);
            if (h != NO_DATA_VALUE)
            {
                VerticalDatum::transform( from, to, lat, lon, h );
            }
        }
    }

    return true;
}

osg::BoundingSphered
FeatureModelGraph::getBoundInWorldCoords(const GeoExtent& extent,
                                         const MapFrame*  mapf ) const
{
    osg::Vec3d center, corner;
    //double z = 0.0;
    GeoExtent workingExtent;

    if ( extent.getSRS()->isEquivalentTo( _usableMapExtent.getSRS() ) )
    {
        workingExtent = extent;
    }
    else
    {
        workingExtent = extent.transform( _usableMapExtent.getSRS() ); // safe.
    }

    workingExtent.getCentroid( center.x(), center.y() );
    
    double centerZ = 0.0;    
    if ( mapf )
    {
        // Use an appropriate resolution for this extents width
        double resolution = workingExtent.width();             
        ElevationQuery query( *mapf );
        query.getElevation( GeoPoint(mapf->getProfile()->getSRS(),center), center.z(), resolution );
        centerZ = center.z();
    }    

    corner.x() = workingExtent.xMin();
    corner.y() = workingExtent.yMin();
    corner.z() = 0;

    if ( _session->getMapInfo().isGeocentric() )
    {
        workingExtent.getSRS()->transformToECEF( center, center );
        workingExtent.getSRS()->transformToECEF( corner, corner );
    }

    return osg::BoundingSphered( center, (center-corner).length() );
}

void
HeightFieldUtils::resolveInvalidHeights(osg::HeightField* grid,
                                        const GeoExtent&  ex,
                                        float             invalidValue,
                                        const Geoid*      geoid)
{
    if ( geoid )
    {
        // need the lat/long extent for geoid queries:
        unsigned numRows = grid->getNumRows();
        unsigned numCols = grid->getNumColumns();
        GeoExtent geodeticExtent = ex.getSRS()->isGeographic() ? ex : ex.transform( ex.getSRS()->getGeographicSRS() );
        double latMin = geodeticExtent.yMin();
        double lonMin = geodeticExtent.xMin();
        double lonInterval = geodeticExtent.width() / (double)(numCols-1);
        double latInterval = geodeticExtent.height() / (double)(numRows-1);

        for( unsigned r=0; r<numRows; ++r )
        {
            double lat = latMin + latInterval*(double)r;
            for( unsigned c=0; c<numCols; ++c )
            {
                double lon = lonMin + lonInterval*(double)c;
                if ( grid->getHeight(c, r) == invalidValue )
                {
                    grid->setHeight( c, r, geoid->getHeight(lat, lon) );
                }
            }
        }
    }
    else
    {
        for(unsigned int i=0; i<grid->getHeightList().size(); i++ )
        {
            if ( grid->getHeightList()[i] == invalidValue )
            {
                grid->getHeightList()[i] = 0.0;
            }
        }
    }
}

GeoImage
GeoImage::crop( const GeoExtent& extent, bool exact, unsigned int width, unsigned int height  ) const
{
    //Check for equivalence
    if ( extent.getSRS()->isEquivalentTo( getSRS() ) )
    {
        //If we want an exact crop or they want to specify the output size of the image, use GDAL
        if (exact || width != 0 || height != 0 )
        {
            OE_DEBUG << "[osgEarth::GeoImage::crop] Performing exact crop" << std::endl;

            //Suggest an output image size
            if (width == 0 || height == 0)
            {
                double xRes = (getExtent().xMax() - getExtent().xMin()) / (double)_image->s();
                double yRes = (getExtent().yMax() - getExtent().yMin()) / (double)_image->t();

                width =  osg::maximum(1u, (unsigned int)((extent.xMax() - extent.xMin()) / xRes));
                height = osg::maximum(1u, (unsigned int)((extent.yMax() - extent.yMin()) / yRes));

                OE_DEBUG << "[osgEarth::GeoImage::crop] Computed output image size " << width << "x" << height << std::endl;
            }

            //Note:  Passing in the current SRS simply forces GDAL to not do any warping
            return reproject( getSRS(), &extent, width, height);
        }
        else
        {
            OE_DEBUG << "[osgEarth::GeoImage::crop] Performing non-exact crop " << std::endl;
            //If an exact crop is not desired, we can use the faster image cropping code that does no resampling.
            double destXMin = extent.xMin();
            double destYMin = extent.yMin();
            double destXMax = extent.xMax();
            double destYMax = extent.yMax();

            osg::Image* new_image = ImageUtils::cropImage(
                _image.get(),
                _extent.xMin(), _extent.yMin(), _extent.xMax(), _extent.yMax(),
                destXMin, destYMin, destXMax, destYMax );

            //The destination extents may be different than the input extents due to not being able to crop along pixel boundaries.
            return new_image?
                GeoImage( new_image, GeoExtent( getSRS(), destXMin, destYMin, destXMax, destYMax ) ) :
                GeoImage::INVALID;
        }
    }
    else
    {
        //TODO: just reproject the image before cropping
        OE_NOTICE << "[osgEarth::GeoImage::crop] Cropping extent does not have equivalent SpatialReference" << std::endl;
        return GeoImage::INVALID;
    }
}

void
MPTerrainEngineNode::invalidateRegion(const GeoExtent& extent,
                                      unsigned         minLevel,
                                      unsigned         maxLevel)
{
    if (_terrainOptions.incrementalUpdate() == true && _liveTiles.valid())
    {
        GeoExtent extentLocal = extent;
        if ( !extent.getSRS()->isEquivalentTo(this->getMap()->getSRS()) )
        {
            extent.transform(this->getMap()->getSRS(), extentLocal);
        }
        
        _liveTiles->setDirty(extentLocal, minLevel, maxLevel);
    }
}

GeoExtent
Profile::clampAndTransformExtent( const GeoExtent& input, bool* out_clamped ) const
{
    if ( out_clamped )
        *out_clamped = false;

    // do the clamping in LAT/LONG.
    const SpatialReference* geo_srs = getSRS()->getGeographicSRS();

    // get the input in lat/long:
    GeoExtent gcs_input =
        input.getSRS()->isGeographic()?
        input :
        input.transform( geo_srs );

    // bail out on a bad transform:
    if ( !gcs_input.isValid() )
        return GeoExtent::INVALID;

    // bail out if the extent's do not intersect at all:
    if ( !gcs_input.intersects(_latlong_extent, false) )
        return GeoExtent::INVALID;

    // clamp it to the profile's extents:
    GeoExtent clamped_gcs_input = GeoExtent(
        gcs_input.getSRS(),
        osg::clampBetween( gcs_input.xMin(), _latlong_extent.xMin(), _latlong_extent.xMax() ),
        osg::clampBetween( gcs_input.yMin(), _latlong_extent.yMin(), _latlong_extent.yMax() ),
        osg::clampBetween( gcs_input.xMax(), _latlong_extent.xMin(), _latlong_extent.xMax() ),
        osg::clampBetween( gcs_input.yMax(), _latlong_extent.yMin(), _latlong_extent.yMax() ) );

    if ( out_clamped )
        *out_clamped = (clamped_gcs_input != gcs_input);

    // finally, transform the clamped extent into this profile's SRS and return it.
    GeoExtent result =
        clamped_gcs_input.getSRS()->isEquivalentTo( this->getSRS() )?
        clamped_gcs_input :
        clamped_gcs_input.transform( this->getSRS() );

    if (result.isValid())
    {
        OE_DEBUG << LC << "clamp&xform: input=" << input.toString() << ", output=" << result.toString() << std::endl;
    }

    return result;
}

bool TileIndex::add( const std::string& filename, const GeoExtent& extent )
{       
    osg::ref_ptr< Polygon > polygon = new Polygon();
    polygon->push_back( osg::Vec3d(extent.bounds().xMin(), extent.bounds().yMin(), 0) );
    polygon->push_back( osg::Vec3d(extent.bounds().xMax(), extent.bounds().yMin(), 0) );
    polygon->push_back( osg::Vec3d(extent.bounds().xMax(), extent.bounds().yMax(), 0) );
    polygon->push_back( osg::Vec3d(extent.bounds().xMin(), extent.bounds().yMax(), 0) );
    polygon->push_back( osg::Vec3d(extent.bounds().xMin(), extent.bounds().yMin(), 0) );
   
    osg::ref_ptr< Feature > feature = new Feature( polygon.get(), extent.getSRS()  );
    feature->set("location", filename );
    
    const SpatialReference* wgs84 = SpatialReference::create("epsg:4326");
    feature->transform( wgs84 );

    return _features->insertFeature( feature.get() );    
    return true;
}

bool
TerrainLayer::mayHaveDataInExtent(const GeoExtent& ex) const
{
    if (!ex.isValid())
    {
        // bad extent; no data
        return false;
    }
    
    const DataExtentList& de = getDataExtents();
    if (de.empty())
    {
        // not enough info, assume yes
        return true;
    }

    // Get extent in local profile:
    GeoExtent localExtent = ex;
    if (getProfile() && getProfile()->getSRS()->isHorizEquivalentTo(ex.getSRS()))
    {
        localExtent = getProfile()->clampAndTransformExtent(ex);
    }

    // Check union:
    if (getDataExtentsUnion().intersects(localExtent))
    {
        // possible yes
        return true;
    }

    // Check each extent in turn:
    for (DataExtentList::const_iterator i = de.begin(); i != de.end(); ++i)
    {
        if (i->intersects(localExtent))
        {
            // possible yes
            return true;
        }
    }

    // definite no.
    return false;
}