← Back to branch summary

~ubuntu-branches/debian/jessie/eso-midas/jessie

« back to all changes in this revision

Viewing changes to prim/general/libsrc/polfil.for

  • Committer: Package Import Robot
  • Author(s): Ole Streicher
  • Date: 2014-04-22 14:44:58 UTC
  • Revision ID: package-import@ubuntu.com-20140422144458-okiwi1assxkkiz39
Tags: upstream-13.09pl1.2+dfsg
Import upstream version 13.09pl1.2+dfsg

Show diffs side-by-side

added added

removed removed

Lines of Context:
prim/general/libsrc/polfil.for
 
1
C @(#)polfil.for        19.1 (ES0-DMD) 02/25/03 14:01:08
 
2
C===========================================================================
 
3
C Copyright (C) 1995 European Southern Observatory (ESO)
 
4
C
 
5
C This program is free software; you can redistribute it and/or 
 
6
C modify it under the terms of the GNU General Public License as 
 
7
C published by the Free Software Foundation; either version 2 of 
 
8
C the License, or (at your option) any later version.
 
9
C
 
10
C This program is distributed in the hope that it will be useful,
 
11
C but WITHOUT ANY WARRANTY; without even the implied warranty of
 
12
C MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 
13
C GNU General Public License for more details.
 
14
C
 
15
C You should have received a copy of the GNU General Public 
 
16
C License along with this program; if not, write to the Free 
 
17
C Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, 
 
18
C MA 02139, USA.
 
19
C
 
20
C Corresponding concerning ESO-MIDAS should be addressed as follows:
 
21
C       Internet e-mail: midas@eso.org
 
22
C       Postal address: European Southern Observatory
 
23
C                       Data Management Division 
 
24
C                       Karl-Schwarzschild-Strasse 2
 
25
C                       D 85748 Garching bei Muenchen 
 
26
C                       GERMANY
 
27
C===========================================================================
 
28
C
 
29
      SUBROUTINE POLFIL(CASE,INIM,TSTIM,IN2IM,
 
30
     +                  OUTIM,NPIX,WINDOW,THRESH,FILL)
 
31
C
 
32
C ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
33
C       
 
34
C.IDENTIFICATION
 
35
C  subroutine POLFIL                    version 1.00    880912
 
36
C  subroutine ISCAN                     version 1.00    880912
 
37
C  subroutine OVRLAP                    version 1.30    881028
 
38
C  subroutine PROGRS                    version 2.20    880919
 
39
C  subroutine MEMDSK                    version 1.00    891129
 
40
C  K. Banse                             ESO - Garching
 
41
C 
 
42
C.KEYWORDS
 
43
C  polygon fill
 
44
C 
 
45
C.PURPOSE
 
46
C  fill pixels inside given polygon (including borderline pixels)
 
47
C 
 
48
C.ALGORITHM
 
49
C  get beginning + end of polygon per line in UP, DOWN
 
50
C  replace all pixels inside with constant or other pixel values
 
51
C       
 
52
C.IN/OUTPUT:
 
53
C call as   POLFIL(CASE,INIM,TSTIM,IN2IM,OUTIM,NPIX,WINDOW,THRESH,FILL)
 
54
C
 
55
C input par:
 
56
C CASE:         char. exp.      = C for constant, I for image
 
57
C INIM:         R*4             input image
 
58
C TSTIM:        R*4             test image
 
59
C IN2IM:        R*4             second input image
 
60
C OUTIM:        R*4             output image
 
61
C NPIX:         I*4 array       NPIX of images
 
62
C WINDOW:       I*4 array       lower left + upper right corner of window
 
63
C                               in images
 
64
C THRESH:       R*4 array       low,high threshold for polygon detection
 
65
C FILL:         R*4             fill value
 
66
C
 
67
C.VERSIONS
 
68
C  1.00         from version 1.00 as of 860707
 
69
C
 
70
C-----------------------------------------------------------------------
 
71
C
 
72
      IMPLICIT NONE
 
73
C
 
74
      CHARACTER*1    CASE
 
75
C
 
76
      INTEGER     WINDOW(4),NPIX(2)
 
77
      INTEGER     NOLINS,NOPIX,OFF,UP,DOWN,SW,KX,JX
 
78
      INTEGER     M,NN,N,KLO,KHI,MLO,MHI,JLO,JHI
 
79
      INTEGER     LOWUP
 
80
C      
 
81
      REAL        INIM(*),TSTIM(*),OUTIM(*),IN2IM(*)
 
82
      REAL        THRESH(2),FILL
 
83
C      
 
84
C  initalize pointer to lower left corner of window
 
85
      OFF = ( (WINDOW(2) - 1) * NPIX(1) ) + WINDOW(1) - 1
 
86
      NOLINS = WINDOW(4) - WINDOW(2) + 1
 
87
      NOPIX = WINDOW(3) - WINDOW(1) + 1
 
88
C
 
89
C  here we go
 
90
C      
 
91
      SW = 0
 
92
      IF (CASE.EQ.'C') THEN
 
93
         DO 1000 N=1,NOLINS
 
94
            IF (SW.NE.0) WRITE(*,12000) (N-1)
 
95
            SW = 0
 
96
            UP = 0
 
97
            DOWN = 0
 
98
            KX = 0
 
99
            JX = 0
 
100
            JLO = 0
 
101
            JHI = 0
 
102
            DO 400, NN=OFF+1,OFF+NOPIX
 
103
               M = NN - 1
 
104
               IF ( (TSTIM(NN).GE.THRESH(1)) .AND. 
 
105
     +              (TSTIM(NN).LE.THRESH(2)) ) THEN
 
106
                  IF (UP.NE.M) THEN
 
107
                     IF (DOWN.EQ.0) THEN
 
108
                        SW = 1 
 
109
                        KLO = LOWUP(TSTIM,THRESH,NPIX,1,NN)
 
110
                        KHI = LOWUP(TSTIM,THRESH,NPIX,2,NN)
 
111
                        KX = NN
 
112
                     ELSE
 
113
                        JX = NN
 
114
                        JLO = LOWUP(TSTIM,THRESH,NPIX,1,NN)
 
115
                        JHI = LOWUP(TSTIM,THRESH,NPIX,2,NN)
 
116
                     ENDIF
 
117
                  ENDIF
 
118
                  UP = NN
 
119
               ELSE
 
120
                  IF (UP.EQ.M) THEN
 
121
                     MLO = LOWUP(TSTIM,THRESH,NPIX,1,M)    !check line below
 
122
                     MHI = LOWUP(TSTIM,THRESH,NPIX,2,M)    !check line above
 
123
                     IF (DOWN.EQ.0) THEN
 
124
                        DOWN = 1 - DOWN
 
125
                        IF (KX.NE.M) THEN
 
126
                           IF ((MLO.EQ.1).AND.(MHI.EQ.1)) THEN
 
127
                              SW = 0
 
128
                              DOWN = 0
 
129
                           ELSE IF ((KLO.EQ.MLO) .AND. 
 
130
     +                              (KHI.EQ.MHI)) THEN
 
131
                              SW = 1 - SW 
 
132
                              DOWN = 0
 
133
                           ENDIF
 
134
                        ELSE IF (MHI.NE.MLO) THEN
 
135
                           IF ((KLO.EQ.MLO) .AND. (KHI.EQ.MHI)) THEN
 
136
                              SW = 1 - SW 
 
137
                              DOWN = 0
 
138
                           ENDIF
 
139
                        ENDIF
 
140
                     ELSE
 
141
                        IF ((JX.NE.M) .OR.(MHI.NE.MLO)) THEN
 
142
                           IF ((JLO.NE.MLO) .OR. (JHI.NE.MHI)) THEN
 
143
                              SW = 1 - SW
 
144
                              DOWN = 0
 
145
                           ENDIF
 
146
                        ELSE
 
147
                           SW = 1 - SW
 
148
                           DOWN = 0
 
149
                        ENDIF
 
150
                     ENDIF
 
151
                  ENDIF
 
152
               ENDIF
 
153
C 
 
154
               IF (SW.NE.0) OUTIM(NN) = FILL
 
155
400         CONTINUE
 
156
            OFF = OFF + NPIX(1)
 
157
1000     CONTINUE
 
158
C      
 
159
C  use second input image
 
160
      ELSE
 
161
         DO 2000 N=1,NOLINS
 
162
            IF (SW.NE.0) WRITE(*,12000) (N-1)
 
163
            SW = 0
 
164
            UP = 0
 
165
            DOWN = 0
 
166
            KX = 0
 
167
            JX = 0
 
168
            JLO = 0
 
169
            JHI = 0
 
170
            DO 1400, NN=OFF+1,OFF+NOPIX
 
171
               M = NN - 1
 
172
               IF ( (TSTIM(NN).GE.THRESH(1)) .AND.
 
173
     +              (TSTIM(NN).LE.THRESH(2)) ) THEN
 
174
                  IF (UP.NE.M) THEN
 
175
                     IF (DOWN.EQ.0) THEN
 
176
                        SW = 1
 
177
                        KLO = LOWUP(TSTIM,THRESH,NPIX,1,NN)
 
178
                        KHI = LOWUP(TSTIM,THRESH,NPIX,2,NN)
 
179
                        KX = NN
 
180
                     ELSE
 
181
                        JX = NN
 
182
                        JLO = LOWUP(TSTIM,THRESH,NPIX,1,NN)
 
183
                        JHI = LOWUP(TSTIM,THRESH,NPIX,2,NN)
 
184
                     ENDIF
 
185
                  ENDIF
 
186
                  UP = NN
 
187
               ELSE
 
188
                  IF (UP.EQ.M) THEN
 
189
                     MLO = LOWUP(TSTIM,THRESH,NPIX,1,M)    !check line below
 
190
                     MHI = LOWUP(TSTIM,THRESH,NPIX,2,M)    !check line above
 
191
                     IF (DOWN.EQ.0) THEN
 
192
                        DOWN = 1 - DOWN
 
193
                        IF (KX.NE.M) THEN
 
194
                           IF ((MLO.EQ.1).AND.(MHI.EQ.1)) THEN
 
195
                              SW = 0
 
196
                              DOWN = 0
 
197
                           ELSE IF ((KLO.EQ.MLO) .AND. 
 
198
     +                              (KHI.EQ.MHI)) THEN
 
199
                              SW = 1 - SW
 
200
                              DOWN = 0
 
201
                           ENDIF
 
202
                        ELSE IF (MHI.NE.MLO) THEN
 
203
                           IF ((KLO.EQ.MLO) .AND. (KHI.EQ.MHI)) THEN
 
204
                              SW = 1 - SW
 
205
                              DOWN = 0
 
206
                           ENDIF
 
207
                        ENDIF
 
208
                     ELSE
 
209
                        IF ((JX.NE.M) .OR.(MHI.NE.MLO)) THEN
 
210
                           IF ((JLO.NE.MLO) .OR. (JHI.NE.MHI)) THEN
 
211
                              SW = 1 - SW
 
212
                              DOWN = 0
 
213
                           ENDIF
 
214
                        ELSE
 
215
                           SW = 1 - SW
 
216
                           DOWN = 0
 
217
                        ENDIF
 
218
                     ENDIF
 
219
                  ENDIF
 
220
               ENDIF
 
221
C
 
222
               IF (SW.NE.0) OUTIM(NN) = IN2IM(NN)
 
223
1400        CONTINUE
 
224
            OFF = OFF + NPIX(1)
 
225
2000     CONTINUE
 
226
C
 
227
      ENDIF
 
228
C 
 
229
C  that's it folks
 
230
      RETURN
 
231
C 
 
232
12000 FORMAT(' line',I5,' has problems...')
 
233
C 
 
234
      END
 
235
 
 
236
      INTEGER FUNCTION LOWUP(A,THR,NPIX,LFLAG,PIXNO)
 
237
C      
 
238
      IMPLICIT NONE
 
239
C      
 
240
      REAL        A(*),THR(2)
 
241
C      
 
242
      INTEGER     PIXNO,NPIX(2),LFLAG
 
243
      INTEGER     N,K
 
244
C      
 
245
      LOWUP = 0
 
246
      N = PIXNO - 1
 
247
      IF (LFLAG.EQ.1) THEN                  !check lower line
 
248
         N = N - NPIX(1)
 
249
         IF (N.LT.1) RETURN
 
250
C 
 
251
         DO 200, K=1,3
 
252
            IF ((A(N).GE.THR(1)).AND.(A(N).LE.THR(2))) GOTO 1000
 
253
            N = N + 1
 
254
200      CONTINUE
 
255
      ELSE
 
256
         N = N + NPIX(1)
 
257
         IF (N.GT.(NPIX(1)*NPIX(2))) RETURN
 
258
C 
 
259
         DO 400, K=1,3
 
260
            IF ((A(N).GE.THR(1)).AND.(A(N).LE.THR(2))) GOTO 1000
 
261
            N = N + 1
 
262
400      CONTINUE
 
263
      ENDIF
 
264
      RETURN
 
265
C      
 
266
1000  LOWUP = 1
 
267
      RETURN
 
268
      END
 
269
 
 
270
      SUBROUTINE OVRLAP(STARTA,STEPA,NPIXA,STARTB,STEPB,NPIXB,
 
271
     +                  BEGIN,END,STAT)
 
272
C
 
273
C++++++++++++++++++++++++++++++++++++++++++++++++++
 
274
C
 
275
C.IDENTIFICATION
 
276
C  subroutine OVRLAP            version 1.30      881028
 
277
C  K. Banse                     ESO - Garching
 
278
C
 
279
C.KEYWORDS
 
280
C  conjunction
 
281
C
 
282
C.PURPOSE
 
283
C  determine overlapping region of two "lines" in real space
 
284
C
 
285
C.ALGORITHM
 
286
C  straight forward
 
287
C
 
288
C.INPUT/OUTPUT
 
289
C  call as   OVRLAP(STARTA,STEPA,NPIXA,STARTB,STEPB,NPIXB,BEGIN,END,STAT)
 
290
C
 
291
C  input par:
 
292
C  STARTA:      R*8     start of "line" A
 
293
C  STEPA:       R*8     stepsize of A
 
294
C  NPIXB:       I*4     no. of steps in A
 
295
C  STARTB:      R*8     start of "line" B
 
296
C  STEPB:       R*8     stepsize of B
 
297
C  NPIXB:       I*4     no. of steps in B
 
298
C
 
299
C  output par:
 
300
C  BEGIN:       R*8     begin of overlapping interval
 
301
C  END:         R*8     end of overlapping interval
 
302
C  STAT:        I*4     return status, should be = 0
 
303
C
 
304
C--------------------------------------------------
 
305
C
 
306
      IMPLICIT NONE
 
307
C
 
308
      INTEGER      NPIXA,NPIXB,STAT
 
309
C
 
310
      DOUBLE PRECISION STARTA,STEPA,STARTB,STEPB
 
311
      DOUBLE PRECISION EB,BEGIN,END
 
312
C
 
313
C  compute position of last pixel
 
314
      STAT = 0
 
315
      END = STARTA + (NPIXA-1)*STEPA
 
316
      EB = STARTB + (NPIXB-1)*STEPB
 
317
C
 
318
C  negative stepsize
 
319
      IF (STEPA.LT.0.D0) THEN
 
320
         IF (STARTA.LT.STARTB) THEN            !low begin
 
321
            BEGIN = STARTA
 
322
         ELSE
 
323
            BEGIN = STARTB
 
324
         ENDIF
 
325
         IF (END.LT.EB) END = EB            !high end
 
326
C
 
327
         IF (END.GT.BEGIN) STAT = 1            !no overlap...!
 
328
C
 
329
C  positive stepsize
 
330
      ELSE
 
331
         IF (STARTA.LT.STARTB) THEN            !high begin
 
332
            BEGIN = STARTB
 
333
         ELSE
 
334
            BEGIN = STARTA
 
335
         ENDIF
 
336
         IF (END.GT.EB) END = EB            !low end
 
337
C
 
338
         IF (END.LT.BEGIN) STAT = 1            !no overlap...!
 
339
      ENDIF
 
340
C
 
341
      RETURN
 
342
      END
 
343
 
 
344
      SUBROUTINE PROGRS(INC,CINC,PRCNT,CHECK)
 
345
C
 
346
C++++++++++++++++++++++++++++++++++++++++++++++++++
 
347
C
 
348
C.IDENTIFICATION
 
349
C  subroutine PROGRS            version 2.00    860609
 
350
C  K. Banse                     ESO - Garching
 
351
C                                       2.20    880919
 
352
C
 
353
C.KEYWORDS
 
354
C  time, progress of work
 
355
C
 
356
C.PURPOSE
 
357
C  display current time + percentage of work already done
 
358
C
 
359
C.ALGORITHM
 
360
C  straight forward 
 
361
C
 
362
C.INPUT/OUTPUT
 
363
C  call as PROGRS(INC,CINC,PRCNT,CHECK)
 
364
C
 
365
C  input par:
 
366
C  INC:         I*4             no. of percent to increase
 
367
C  CINC:        I*4             increment for check
 
368
C  
 
369
C  in/output par:
 
370
C  PRCNT:       I*4             percent of processing done
 
371
C  CHECK:       I*4             check value
 
372
C 
 
373
C.VERSIONS
 
374
C  2.10         do not write to logfile - only to terminal...
 
375
C  2.20         move to FORTRAN 77 + use system independent calls
 
376
C 
 
377
C--------------------------------------------------
 
378
C
 
379
      IMPLICIT NONE
 
380
C      
 
381
      INTEGER     INC,CINC,PRCNT,CHECK
 
382
      INTEGER     N,IAV,KLOG,NLOG,STAT
 
383
      INTEGER     UNIT(1),NULLO
 
384
C      
 
385
      CHARACTER   TIME*40,PBUF*60
 
386
C      
 
387
      IF (PRCNT.GT.100) RETURN
 
388
C
 
389
C  get current time, convert percentage to ASCII + display all that
 
390
      TIME = ' '
 
391
      CALL GENTIM(TIME)
 
392
      DO 100 N=40,2,-1
 
393
         IF (TIME(N:N).NE.' ') THEN
 
394
            IAV = N
 
395
            GOTO 200
 
396
         ENDIF
 
397
100   CONTINUE
 
398
      IAV = 1
 
399
C 
 
400
200   PBUF(1:) = ' '
 
401
      WRITE(PBUF(1:),10000) TIME(1:IAV),PRCNT
 
402
C      
 
403
C  set log flag to 0, so we do not write into the log file...
 
404
      CALL STKRDI('LOG',1,1,IAV,KLOG,UNIT,NULLO,STAT)
 
405
      NLOG = 0
 
406
      CALL STKWRI('LOG',NLOG,1,1,UNIT,STAT)
 
407
      CALL STTPUT(PBUF,STAT)
 
408
      CALL STKWRI('LOG',KLOG,1,1,UNIT,STAT)
 
409
C      
 
410
C  update
 
411
      PRCNT = PRCNT + INC
 
412
      CHECK = CHECK + CINC
 
413
      RETURN
 
414
C      
 
415
10000 FORMAT(A,I4,'% done ... ')
 
416
C      
 
417
      END
 
418
 
 
419
      SUBROUTINE MEMDSK(IMNO,FLAG,A,NPIX,STA,WSIZE,STAT)
 
420
C
 
421
C++++++++++++++++++++++++++++++++++++++++++++++++++
 
422
C
 
423
C.IDENTIFICATION
 
424
C  subroutine MEMDSK            version 1.00      891129
 
425
C  K. Banse                     ESO - Garching
 
426
C
 
427
C.KEYWORDS
 
428
C  memory, diskfile
 
429
C
 
430
C.PURPOSE
 
431
C  copy data stored in virtual memory to right place in disk file
 
432
C
 
433
C.ALGORITHM
 
434
C  use STFPUT interface
 
435
C
 
436
C.INPUT/OUTPUT
 
437
C  call as   MEMDSK(IMNO,FLAG,A,NPIX,STA,WSIZE,STAT)
 
438
C
 
439
C  input par:
 
440
C  IMNO:        Integer         file no.
 
441
C  FLAG:        Integer         = 0, buffer A holds full frame
 
442
C                               = 1, buffer A holds data from window start on
 
443
C                               = 2, buffer A holds only the data of the window
 
444
C  A:           Real            buffer with data
 
445
C  NPIX:        Integer         no. of pixels in A
 
446
C  STA:         Integer         start of window in total frame
 
447
C  WSIZE:       Integer         no. of pixels in window
 
448
C
 
449
C  output par:
 
450
C  STAT:        Integer         return status, should be = 0
 
451
C
 
452
C--------------------------------------------------
 
453
C
 
454
      IMPLICIT NONE
 
455
C 
 
456
      INTEGER     IMNO,FLAG,NPIX(*),STA(*),WSIZE(*),STAT
 
457
      INTEGER     FELEM,SIZE,N,IOFF
 
458
C 
 
459
      REAL        A(*)
 
460
C 
 
461
      FELEM = ((STA(2)-1) * NPIX(1)) + STA(1)
 
462
      IF (FLAG .EQ. 0) THEN
 
463
         SIZE = WSIZE(2) * NPIX(1)              !WSIZE(2) lines...
 
464
         CALL STFPUT(IMNO,FELEM,SIZE,A(FELEM),STAT)
 
465
      ELSE IF (FLAG .EQ. 1) THEN
 
466
         SIZE = WSIZE(2) * NPIX(1)              !WSIZE(2) lines...
 
467
         CALL STFPUT(IMNO,FELEM,SIZE,A(1),STAT)
 
468
      ELSE IF (FLAG .EQ. 2) THEN
 
469
         IOFF = 1
 
470
         DO 1000 N=1,WSIZE(2)                   !again WSIZE(2) lines ...
 
471
            CALL STFPUT(IMNO,FELEM,WSIZE(1),A(IOFF),STAT)
 
472
            IOFF = IOFF + WSIZE(1)
 
473
            FELEM = FELEM + NPIX(1)
 
474
1000     CONTINUE
 
475
      ENDIF
 
476
C 
 
477
      RETURN
 
478
      END
 
479