Skip to content
Snippets Groups Projects
profilevisualization.py 84.4 KiB
Newer Older
  • Learn to ignore specific revisions
  • Benjamin Jakimow's avatar
    Benjamin Jakimow committed
            Loads a temporal profile for a single or multiple geometries.
            :param spatialPoints: SpatialPoint | [list-of-SpatialPoints]
            """
    
            if not isinstance(self.plotSettingsModel2D, PlotSettingsModel2D):
    
            #if not self.pixelLoader.isReadyToLoad():
            #    return False
    
    
            assert isinstance(self.TS, TimeSeries)
    
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
            #Get or create the TimeSeriesProfiles which will store the loaded values
    
            tasks = []
            TPs = []
            theGeometries = []
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
    
            # Define a which (new) bands need to be loaded for each sensor
    
            if LUT_bandIndices is None:
                LUT_bandIndices = dict()
                for sensor in self.TS.Sensors:
    
                    LUT_bandIndices[sensor] = self.plotSettingsModel2D.requiredBandsIndices(sensor)
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
    
    
            assert isinstance(LUT_bandIndices, dict)
            for sensor in self.TS.Sensors:
                assert sensor in LUT_bandIndices.keys()
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
    
            #update new / existing points
            if isinstance(spatialPoints, SpatialPoint):
                spatialPoints = [spatialPoints]
    
            for spatialPoint in spatialPoints:
                assert isinstance(spatialPoint, SpatialPoint)
                TP = self.tpCollection.fromSpatialPoint(spatialPoint)
                if not isinstance(TP, TemporalProfile):
                    TP = TemporalProfile(self.TS, spatialPoint)
    
                    self.tpCollection.insertTemporalProfiles(TP, i=0)
    
    
                    if len(self.tpCollection) == 1:
                        if len(self.plotSettingsModel2D) == 0:
    
                            self.createNewPlotStyle2D()
    
    
                        if len(self.plotSettingsModel3D) == 0:
                            #todo: individual 3D style
                            pass
    
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                TPs.append(TP)
                theGeometries.append(TP.mCoordinate)
    
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
            TP_ids = [TP.mID for TP in TPs]
            #each TSD is a Task
            #a Task defines which bands are to be loaded
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                #do not load from invisible TSDs
                if not tsd.isVisible():
                    continue
    
                #which bands do we need to load?
                requiredIndices = set(LUT_bandIndices[tsd.sensor])
                if len(requiredIndices) == 0:
                    continue
                else:
                    s = ""
    
                missingIndices = set()
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                for TP in TPs:
                    assert isinstance(TP, TemporalProfile)
    
                    need2load = TP.missingBandIndices(tsd, requiredIndices=requiredIndices)
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                    missingIndices = missingIndices.union(need2load)
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                missingIndices = sorted(list(missingIndices))
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                if len(missingIndices) > 0:
                    task = PixelLoaderTask(tsd.pathImg, theGeometries,
                                           bandIndices=missingIndices,
                                           temporalProfileIDs=TP_ids)
                    tasks.append(task)
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
            if len(tasks) > 0:
                aGoodDefault = 2 if len(self.TS) > 25 else 1
    
    
                #self.pixelLoader.setNumberOfProcesses(SETTINGS.value('profileloader_threads', aGoodDefault))
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                if DEBUG:
                    print('Start loading for {} geometries from {} sources...'.format(
                        len(theGeometries), len(tasks)
                    ))
                self.pixelLoader.startLoading(tasks)
    
            else:
                if DEBUG:
                    print('Data for geometries already loaded')
    
        def setVisibility(self, sensorPlotStyle):
    
            assert isinstance(sensorPlotStyle, TemporalProfile2DPlotStyle)
    
            self.setVisibility2D(sensorPlotStyle)
    
        def setVisibility2D(self, sensorPlotStyle):
    
        def addData(self, sensorView = None):
    
                for sv in self.plotSettingsModel2D.items:
    
                assert isinstance(sensorView, TemporalProfile2DPlotStyle)
    
        @QtCore.pyqtSlot()
    
            for plotSetting in self.plotSettingsModel2D:
                assert isinstance(plotSetting, TemporalProfile2DPlotStyle)
    
                tp = plotSetting.temporalProfile()
                for pdi in plotSetting.mPlotItems:
                    assert isinstance(pdi, TemporalProfilePlotDataItem)
                    pdi.updateDataAndStyle()
                if isinstance(tp, TemporalProfile) and plotSetting.temporalProfile().updated():
    
            for i in self.plot2D.getPlotItem().dataItems:
                i.updateItems()
    
    
    
            notInit = [0, 1] == self.plot2D.getPlotItem().getAxis('bottom').range
            if notInit:
    
                for plotSetting in self.plotSettingsModel2D:
                    assert isinstance(plotSetting, TemporalProfile2DPlotStyle)
    
                    for pdi in plotSetting.mPlotItems:
                        assert isinstance(pdi, TemporalProfilePlotDataItem)
    
                        if pdi.xData.ndim == 0 or pdi.xData.shape[0] == 0:
                            continue
    
                        if x0 is None:
                            x0 = pdi.xData.min()
                            x1 = pdi.xData.max()
                        else:
                            x0 = min(pdi.xData.min(), x0)
                            x1 = max(pdi.xData.max(), x1)
    
                if x0 is not None:
                    self.plot2D.getPlotItem().setXRange(x0, x1)
    
                    #self.plot2D.xAxisInitialized = True
    
    
        @QtCore.pyqtSlot()
        def updatePlot3D(self):
            if OPENGL_AVAILABLE:
                from pyqtgraph.opengl import GLViewWidget
    
                from pyqtgraph.opengl.GLGraphicsItem import GLGraphicsItem
    
                import pyqtgraph.opengl as gl
                assert isinstance(self.plot3D, GLViewWidget)
    
    
    
    
                # 1. ensure that data from all bands will be loaded
                LUT_bandIndices = dict()
                for sensor in self.TS.sensors():
                    assert isinstance(sensor, SensorInstrument)
                    LUT_bandIndices[sensor] = list(range(sensor.nb))
    
                coordinates = []
                for plotStyle3D in self.plotSettingsModel3D:
                    assert isinstance(plotStyle3D, TemporalProfile3DPlotStyle)
                    tp = plotStyle3D.temporalProfile()
                    if isinstance(tp, TemporalProfile):
                        coordinates.append(tp.mCoordinate)
                if len(coordinates) > 0:
                    self.loadCoordinate(coordinates, LUT_bandIndices=LUT_bandIndices)
    
                # 2. remove old plot items
    
                self.ui.plotWidget3D.clearItems()
    
    
                # 3 add new plot items
                for plotStyle3D in self.plotSettingsModel3D:
                    assert isinstance(plotStyle3D, TemporalProfile3DPlotStyle)
                    gli = plotStyle3D.createPlotItem(None)
    
                    if isinstance(gli, GLGraphicsItem):
    
                        self.ui.plotWidget3D.addItem(gli)
    
    
                # w.setBackgroundColor(QColor('black'))
                # w.setCameraPosition(pos=(0.0, 0.0, 0.0), distance=1.)
    
                #self.plot3D.addItem(self.ui.plotWidget3D.glGridItem)
    
                self.plot3D.updateDataRanges()
    
                self.plot3D.update()
    
                self.plot3D.zoomToFull()
    
    
        @QtCore.pyqtSlot()
        def updatePlot2D(self):
    
            if isinstance(self.plotSettingsModel2D, PlotSettingsModel2D):
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                    print('Update plot...')
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                pi = self.plot2D.getPlotItem()
                piDataItems = pi.listDataItems()
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                locations = set()
    
                for plotSetting in self.plotSettingsModel2D:
                    assert isinstance(plotSetting, TemporalProfile2DPlotStyle)
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                    locations.add(plotSetting.temporalProfile().mCoordinate)
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                    for pdi in plotSetting.mPlotItems:
                        assert isinstance(pdi, TemporalProfilePlotDataItem)
    
                        pdi.updateDataAndStyle()
    
    
                #for i in pi.dataItems:
                #    i.updateItems()
    
                #self.plot2D.update()
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                #2. load pixel data
                self.loadCoordinate(list(locations))
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
                # https://github.com/pyqtgraph/pyqtgraph/blob/5195d9dd6308caee87e043e859e7e553b9887453/examples/customPlot.py
                return
    
    def examplePixelLoader():
    
    
        # prepare QGIS environment
        if sys.platform == 'darwin':
            PATH_QGS = r'/Applications/QGIS.app/Contents/MacOS'
            os.environ['GDAL_DATA'] = r'/usr/local/Cellar/gdal/1.11.3_1/share'
        else:
            # assume OSGeo4W startup
            PATH_QGS = os.environ['QGIS_PREFIX_PATH']
        assert os.path.exists(PATH_QGS)
    
        qgsApp = QgsApplication([], True)
        QApplication.addLibraryPath(r'/Applications/QGIS.app/Contents/PlugIns')
        QApplication.addLibraryPath(r'/Applications/QGIS.app/Contents/PlugIns/qgis')
        qgsApp.setPrefixPath(PATH_QGS, True)
        qgsApp.initQgis()
    
    
        gb = QGroupBox()
        gb.setTitle('Sandbox')
    
        PL = PixelLoader()
    
        PL.setNumberOfThreads(2)
    
    
        if False:
            files = ['observationcloud/testdata/2014-07-26_LC82270652014207LGN00_BOA.bsq',
                     'observationcloud/testdata/2014-08-03_LE72270652014215CUB00_BOA.bsq'
                     ]
        else:
    
            from timeseriesviewer.utils import file_search
    
            searchDir = r'H:\LandsatData\Landsat_NovoProgresso'
            files = file_search(searchDir, '*227065*band4.img', recursive=True)
            #files = files[0:3]
    
        lyr = QgsRasterLayer(files[0])
        coord = lyr.extent().center()
        crs = lyr.crs()
    
        l = QVBoxLayout()
    
        btnStart = QPushButton()
        btnStop = QPushButton()
        prog = QProgressBar()
        tboxResults = QPlainTextEdit()
        tboxResults.setMaximumHeight(300)
        tboxThreads = QPlainTextEdit()
        tboxThreads.setMaximumHeight(200)
        label = QLabel()
        label.setText('Progress')
    
        def showProgress(n,m,md):
            prog.setMinimum(0)
            prog.setMaximum(m)
            prog.setValue(n)
    
            info = []
            for k, v in md.items():
                info.append('{} = {}'.format(k,str(v)))
            tboxResults.setPlainText('\n'.join(info))
            #tboxThreads.setPlainText(PL.threadInfo())
            qgsApp.processEvents()
    
        PL.sigPixelLoaded.connect(showProgress)
        btnStart.setText('Start loading')
        btnStart.clicked.connect(lambda : PL.startLoading(files, coord, crs))
        btnStop.setText('Cancel')
        btnStop.clicked.connect(lambda: PL.cancelLoading())
        lh = QHBoxLayout()
        lh.addWidget(btnStart)
        lh.addWidget(btnStop)
        l.addLayout(lh)
        l.addWidget(prog)
        l.addWidget(tboxThreads)
        l.addWidget(tboxResults)
    
        gb.setLayout(l)
        gb.show()
        #rs.setBackgroundStyle('background: qlineargradient(x1:0, y1:0, x2:0, y2:1, stop:0 #222, stop:1 #333);')
        #rs.handle.setStyleSheet('background: qlineargradient(x1:0, y1:0, x2:0, y2:1, stop:0 #282, stop:1 #393);')
        qgsApp.exec_()
        qgsApp.exitQgis()
    
    if __name__ == '__main__':
        import site, sys
    
        from timeseriesviewer import utils
        qgsApp = utils.initQgisApplication()
    
        DEBUG = False
    
        if False: #the ultimative test for floating point division correctness, at least on a DOY-level
            date1 = np.datetime64('1960-12-31','D')
            assert date1 == num2date(date2num(date1))
            #1960 - 12 - 31
            for year in  range(1960, 2057):
                for doy in range(1, daysPerYear(year)+1):
                    dt = datetime.timedelta(days=doy - 1)
                    date1 = np.datetime64('{}-01-01'.format(year)) + np.timedelta64(doy-1,'D')
                    date2 = datetime.date(year=year, month=1, day=1) + datetime.timedelta(days=doy-1)
    
                    assert date1 == num2date(date2num(date1), dt64=True), 'date1: {}'.format(date1)
                    assert date2 == num2date(date2num(date2), dt64=False), 'date2: {}'.format(date1)
    
        ui = ProfileViewDockUI()
        ui.show()
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
    
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
        if True:
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
            STVis = SpectralTemporalVisualization(ui)
            STVis.setTimeSeries(TS)
    
            import example.Images
            from timeseriesviewer import file_search
            files = file_search(os.path.dirname(example.Images.__file__), '*.tif')
    
            TS.addFiles([files[0]])
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
            ext = TS.getMaxSpatialExtent()
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
            cp1 = SpatialPoint(ext.crs(),ext.center())
    
            cpND = SpatialPoint(ext.crs(), 681151.214,-752388.476)
    
            cp2 = SpatialPoint(ext.crs(), ext.center())
            cp3 = SpatialPoint(ext.crs(), ext.center().x()+500, ext.center().y()+250)
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
    
    
            STVis.loadCoordinate(cp2)
            STVis.loadCoordinate(cp3)
    
            STVis.createNewPlotStyle2D()
    
            if False:
                for tp in STVis.tpCollection:
                    assert isinstance(tp, TemporalProfile)
                    tp.plot()
    
            STVis.tpCollection.removeTemporalProfiles(STVis.tpCollection[-1])
            STVis.createNewPlotStyle3D()
    
            STVis.ui.listWidget.setCurrentRow(1)
    
    Benjamin Jakimow's avatar
    Benjamin Jakimow committed
        qgsApp.exec_()
        qgsApp.exitQgis()