← Back to branch summary

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

« back to all changes in this revision

Viewing changes to contrib/romafot/src/rfotselect.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:
contrib/romafot/src/rfotselect.for
 
1
C @(#)rfotselect.for    19.1 (ES0-DMD) 02/25/03 13:30:15
 
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
       PROGRAM SELECT
 
30
C++
 
31
C.IDENTIFICATION: RFOTSELECT
 
32
C.PURPOSE: Select objects from a displayed frame and stored them in an
 
33
C          intermediate table
 
34
C.AUTHOR:  R. Buonanno, G. Buscema, C. Corsi, I. Ferraro, G. Iannicola
 
35
C          Osservatorio Astronomico di Roma
 
36
C.VERSION: 870417  First running version at ESO (outside MIDAS)   R. Buonanno
 
37
C.         881031  New version                                    R. Buonanno
 
38
C.         890220  Major changes for portable MIDAS               R.H. Warmels
 
39
C.                 load of frame done by LOAD/IMAGE
 
40
C.                 cursor interaction with GETCURS
 
41
C----
 
42
      IMPLICIT         NONE
 
43
      INCLUDE          'MID_REL_INCL:RFOTDECL.INC/NOLIST'
 
44
C
 
45
C     INTEGER          AREA(4)
 
46
      INTEGER          CCOUNT,COOFF
 
47
      INTEGER          EC, ED, EL
 
48
      INTEGER          IWND(2)    
 
49
      INTEGER          IVX(9)
 
50
      INTEGER*8        IPNTR
 
51
      INTEGER          IMF
 
52
      INTEGER          I, IT, INX
 
53
      INTEGER          IAC, ICOM, IROW, IDUM
 
54
      INTEGER          ISTAT,ISTAT1,ISTAT2
 
55
      INTEGER          J
 
56
      INTEGER          KUN, KNUL, KL
 
57
      INTEGER          KDX, KDY
 
58
      INTEGER          LX1, LY1
 
59
      INTEGER          LX, LY
 
60
      INTEGER          MADRID(1)
 
61
      INTEGER          NPIX(3)
 
62
      INTEGER          NC, NLIV
 
63
      INTEGER          NCP,NHL, NCH
 
64
      INTEGER          NAXIS, NPL,NL
 
65
      INTEGER          XY1,XY2
 
66
C
 
67
      REAL             AL
 
68
      REAL             BB, B1
 
69
      REAL             DIMX, DIMY
 
70
      REAL             FPIX1(2),FPIX2(2)
 
71
      REAL             FG, FO1
 
72
      REAL             H1
 
73
      REAL             PAR(4)
 
74
      REAL             PARS(4)
 
75
      REAL             RINPUT(4)
 
76
      REAL             RX(100),RY(100)
 
77
      REAL             RINTD1, RINTD2, RINTD3, RINTD4
 
78
      REAL             RINTD5, RINTD6, RINTD7
 
79
      REAL             RME
 
80
      REAL             SQR, SQM
 
81
      REAL             SL
 
82
      REAL             X1, Y1
 
83
      REAL             VET(100)
 
84
      REAL             VALUE1,VALUE2
 
85
      REAL             VME
 
86
      REAL             UU, U1
 
87
      REAL             VV
 
88
      REAL             WC01(2),WC02(2)
 
89
C
 
90
      DOUBLE PRECISION BEGIN(3),STEP(3),DDUM
 
91
C
 
92
      LOGICAL          NUL
 
93
     
 
94
      CHARACTER*60     FRAME
 
95
      CHARACTER*72     CUNIT
 
96
      CHARACTER*80     IDENT
 
97
      CHARACTER*60     INTFIL
 
98
      CHARACTER*3      ACTION
 
99
      CHARACTER*80     STRING
 
100
      CHARACTER*60     CINPUT
 
101
C
 
102
C *** this part into the common blocks
 
103
      INTEGER          NCINT, NRINT, NSINT, NACINT, NARINT
 
104
      INTEGER          TIDINT
 
105
      INTEGER          NGRP, NOBJ, NINT
 
106
C 
 
107
      INCLUDE          'MID_INCLUDE:IDIMEM.INC'
 
108
      INCLUDE          'MID_INCLUDE:IDIDEV.INC'
 
109
      INCLUDE          'MID_INCLUDE:ST_DEF.INC'
 
110
      COMMON           /VMR/MADRID
 
111
      INCLUDE          'MID_INCLUDE:ST_DAT.INC'
 
112
C
 
113
      DATA             IVX/1,2,3,4,5,6,7,8,9/
 
114
      DATA             NCP/1/
 
115
      DATA             NHL/0/
 
116
C
 
117
 9001 FORMAT('*** INFO: Intermediate table contains', I5,
 
118
     2       ' groups with ',I5,' components')
 
119
 9004 FORMAT('          table will be appended')
 
120
 9003 FORMAT('Enter identification, mag., col1 and col2 [',I6,
 
121
     2       ',0.0,0.0,0.0]: ')
 
122
C
 
123
C *** start the code
 
124
      CALL STSPRO('SELECT')
 
125
C
 
126
C *** input parameters
 
127
      CALL STKRDC('IN_A',1,1,60,IAC,FRAME,KUN,KNUL,ISTAT)      ! name of image
 
128
      CALL STIGET(FRAME,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE,
 
129
     2             2,NAXIS,NPIX,BEGIN,STEP,IDENT,CUNIT,IPNTR,IMF,ISTAT)
 
130
      NPL = NPIX(1)
 
131
      NL  = NPIX(2)
 
132
C
 
133
C *** get the intemediate file name
 
134
      CALL STKRDC('OUT_A',1,1,60,IAC,INTFIL,KUN,KNUL,ISTAT)   ! intermediate file
 
135
C
 
136
C *** open or create the intermediate table
 
137
      CALL STECNT('GET',EC,ED,EL)
 
138
      CALL STECNT('PUT',1,0,0)
 
139
      CALL TBTOPN(INTFIL,F_IO_MODE,TIDINT,ISTAT)
 
140
      IF (ISTAT.NE.0) THEN
 
141
         STRING = '*** INFO: Intermediate table not present; '//
 
142
     2            'will create a new one ...'
 
143
         IDNGRP = 0
 
144
         NGRP   = 0
 
145
         NOBJ   = 0
 
146
         IROW   = 0
 
147
         CALL STTPUT(STRING,ISTAT)
 
148
         CALL INTINI(INTFIL,TIDINT)
 
149
 
 
150
      ELSE
 
151
         CALL TBIGET(TIDINT,NCINT,NRINT,NSINT,NACINT,NARINT,ISTAT)
 
152
         IF (ISTAT.NE.0) THEN
 
153
            STRING = '*** FATAL: Problems with getting info for '//
 
154
     2               'intermediate table ...'
 
155
           CALL STTPUT(STRING,ISTAT)
 
156
           CALL STSEPI
 
157
         ENDIF
 
158
C
 
159
         IF (NRINT.EQ.0) THEN
 
160
            STRING = '*** INFO: No points in the intermediate table'
 
161
            CALL STTPUT(STRING,ISTAT)
 
162
         ENDIF
 
163
C
 
164
C *** get the intermediate table parameters
 
165
         CALL INTDRD(TIDINT,NGRP,NOBJ,NINT,RINTD1,RINTD2,RINTD3,
 
166
     2                      RINTD4,RINTD5,RINTD6,RINTD7)
 
167
         IROW = NRINT                                   ! number of rows written
 
168
         WRITE(STRING,9001) NGRP,NOBJ
 
169
         CALL STTPUT(STRING,ISTAT)
 
170
         WRITE(STRING,9004)
 
171
         CALL STTPUT(STRING,ISTAT)
 
172
         CALL TBERDI(TIDINT,NRINT,1,IDNGRP,NUL,ISTAT)
 
173
      ENDIF
 
174
      CALL STECNT('PUT',EC,ED,EL)
 
175
C
 
176
C *** open the display device
 
177
      CALL DTOPEN(1,ISTAT)
 
178
C
 
179
C *** get the area
 
180
C     CALL STKRDC('INPUTC',1,1,40,IAC,STRING,KUN,KNUL,ISTAT)        ! user area
 
181
C     CALL EXTRCO(STRING,NAXIS,NPIX,BEGIN,STEP,NDUM,AREA(1),
 
182
C   2                   AREA(3),ISTAT)
 
183
C     IX0  = AREA(1)
 
184
C     IY0  = AREA(2)
 
185
C     NX   = AREA(3) - AREA(1) + 1
 
186
C     NY   = AREA(4) - AREA(2) + 1
 
187
C
 
188
C     IF (IX0.LT.1) IX0 = 1
 
189
C     IF (IY0.LT.1) IY0 = 1
 
190
C     IF ((IX0+NX-1).GT.NPL) NX = NPL-IX0+1
 
191
C     IF ((IY0+NY-1).GT.NL)  NY = NL-IY0+1
 
192
C     WRITE (STRING,9001) AREA(1),AREA(2),AREA(3),AREA(4)
 
193
C     CALL STTPUT(STRING,ISTAT)
 
194
C
 
195
C *** get the window size
 
196
      CALL STKRDI('INPUTI',1,2,IAC,IWND,KUN,KNUL,ISTAT)
 
197
      DIMX = FLOAT(IWND(1))
 
198
      DIMY = FLOAT(IWND(2))
 
199
      IF (DIMX.GT.55) THEN
 
200
         DIMX = 55
 
201
      ENDIF
 
202
      IF (DIMY.GT.55) THEN
 
203
         DIMY = 55
 
204
      ENDIF
 
205
      KDX = DIMX/2+1
 
206
      KDY = DIMY/2+1
 
207
C
 
208
C *** infinite loop for cursor interaction
 
209
      COOFF     = 0
 
210
      CCOUNT    = 0
 
211
      ACTION    = 'YYN'
 
212
      FRAME     = ' '
 
213
      FLGCMP(1) = 1
 
214
      CALL STTPUT('*** INFO: Use cursor control panel to '//
 
215
     2                       'select object',ISTAT)
 
216
C
 
217
 1000 CONTINUE
 
218
      IDNGRP = IDNGRP + 1
 
219
 1001 CONTINUE
 
220
      CALL GETCUR(ACTION,FRAME,XY1,FPIX1,WC01,VALUE1,ISTAT1,
 
221
     2                         XY2,FPIX2,WC02,VALUE2,ISTAT2)
 
222
      IF (ISTAT1.EQ.0 .AND. ISTAT2. EQ. 0) THEN
 
223
         IF (CCOUNT.EQ.0) THEN
 
224
            IF (COOFF.EQ.1) THEN
 
225
               CALL DTCLOS(QDSPNO)
 
226
               CALL STSEPI
 
227
            ELSE
 
228
               FRAME     = ' '
 
229
               CALL STTPUT('*** WARNING: Switch cursor on ...'//
 
230
     2                     ' next time we exit',ISTAT)
 
231
            ENDIF
 
232
            COOFF = 1
 
233
            GO TO 1001
 
234
 
 
235
         ELSE
 
236
            NINT   = 0
 
237
            RINTD1 = 15000.0
 
238
            RINTD2 = 1.
 
239
            RINTD3 = 3.
 
240
            RINTD4 = 4.
 
241
            RINTD5 = 0.0
 
242
            RINTD6 = 0.0
 
243
            RINTD7 = 0.0
 
244
            CALL INTDWR(TIDINT,NGRP,NOBJ,NINT,RINTD1,RINTD2,RINTD3,
 
245
     2                         RINTD4,RINTD5,RINTD6,RINTD7)
 
246
            CALL TBSINI(TIDINT,ISTAT)
 
247
            CALL TBTCLO(TIDINT,ISTAT)
 
248
            CALL DTCLOS(QDSPNO)
 
249
            CALL STSEPI
 
250
         ENDIF
 
251
 
 
252
      ELSE
 
253
         CCOUNT = 1
 
254
         LX     = INT(FPIX1(1))
 
255
         LY     = INT(FPIX1(2))
 
256
C
 
257
C  *** fit data
 
258
         AL     = 4*ALOG(2.)
 
259
         NC     = 0
 
260
         DO 110 J = 1,9
 
261
            RX(J) = 0
 
262
            RY(J) = 0
 
263
  110    CONTINUE
 
264
C
 
265
         H1 = -10.E15
 
266
         DO 120 I = LY-4,LY+4
 
267
            NC  = NC+1
 
268
            INX = LX-4
 
269
            CALL REALIN(NPL,NL,I,INX,9,MADRID(IPNTR),VET)
 
270
            RME = 0
 
271
C
 
272
            DO 121 J=1,9
 
273
               H1  = AMAX1(H1,VET(J))
 
274
               RME = RME+VET(J)
 
275
               RX(J) = RX(J)+VET(J)
 
276
  121       CONTINUE
 
277
            RY(NC) = RME/9.
 
278
  120    CONTINUE
 
279
C
 
280
         PAR(1) = 0.
 
281
         PAR(4) = 10.E15
 
282
         DO 130 J = 1,9
 
283
            RX(J) = RX(J)/9.
 
284
            IF (RX(J).GT.PAR(1)) THEN
 
285
               PAR(1) = RX(J)
 
286
               PAR(2)=J
 
287
            END IF
 
288
            IF (RX(J).LT.PAR(4)) PAR(4)=RY(J)
 
289
  130    CONTINUE
 
290
C
 
291
         PAR(1) = PAR(1)-PAR(4)
 
292
         KL     = PAR(2)
 
293
         VME    = (RX(KL-1)+RX(KL+1))/2.
 
294
         IF (VME.GE.PAR(1)) VME=PAR(1)/2.
 
295
         PAR(3) = SQRT(AL/ALOG(PAR(1)/VME))
 
296
         SQR    = 10.**15
 
297
         SL     = SQR
 
298
         IT     = 0
 
299
         NLIV   = 9
 
300
C
 
301
         DO 140 J = 1,4
 
302
            PARS(J)=PAR(J)
 
303
  140    CONTINUE
 
304
 
 
305
 2001    CONTINUE
 
306
            IF (SL.LE.0.0001.OR.IT.GT.20) GO TO 2000
 
307
            IT = IT+1
 
308
            CALL MONO4(IVX,RX,NLIV,PAR,.7)
 
309
            IF (PAR(1).LT..1 .OR. ABS(PAR(2)-PARS(2)).GT.4 .OR.
 
310
     2         PAR(3).LE.0. .OR. PAR(3).GT.20.) THEN
 
311
               DO 150 J=1,4
 
312
                  PAR(J) = PARS(J)
 
313
  150          CONTINUE
 
314
               GO TO 2000
 
315
            ENDIF
 
316
C
 
317
            SQM = 0.
 
318
            DO 160  I = 1,NLIV
 
319
               FG  = PAR(1)*EXP(-4*ALOG(2.)*((I-PAR(2))/PAR(3))**2)
 
320
               SQM = SQM+(RX(I)-FG-PAR(4))**2
 
321
  160       CONTINUE
 
322
            SQM    = SQRT(SQM/NLIV)
 
323
            SL     = ABS(SQR-SQM)/SQR
 
324
            SQR    = SQM
 
325
           GO TO 2001
 
326
 2000    CONTINUE
 
327
 
 
328
         X1     = PAR(2)+LX-5
 
329
         FO1    = PAR(4)
 
330
C 
 
331
         PAR(1) = 0
 
332
         PAR(4) = 10.E15
 
333
         DO 170  J=1,9
 
334
            IF (RY(J).GT.PAR(1)) THEN
 
335
               PAR(1) = RY(J)
 
336
               PAR(2) = J
 
337
            END IF
 
338
            IF (RY(J).LT.PAR(4)) THEN
 
339
               PAR(4) = RY(J)
 
340
            ENDIF
 
341
  170    CONTINUE
 
342
C
 
343
         PAR(1) = PAR(1)-PAR(4)
 
344
         KL     = PAR(2)
 
345
         VME    = (RY(KL-1)+RY(KL+1))/2.
 
346
         IF (VME.GE.PAR(1)) VME = PAR(1)/2.
 
347
         PAR(3) = SQRT(AL/ALOG(PAR(1)/VME))
 
348
         SQR    = 10.**15
 
349
         SL     = SQR
 
350
         IT     = 0
 
351
         NLIV   = 9
 
352
C
 
353
         DO 180 J = 1,4
 
354
            PARS(J) = PAR(J)
 
355
  180    CONTINUE
 
356
C
 
357
 3001    CONTINUE
 
358
            IF (SL.LE.0.0001.OR.IT.GT.20) GO TO 3000
 
359
            IT = IT+1
 
360
            CALL MONO4(IVX,RY,NLIV,PAR,.7)
 
361
            IF (PAR(1).LT..1 .OR. ABS(PAR(2)-PARS(2)).GT.4 .OR. 
 
362
     2          PAR(3).LE.0. .OR. PAR(3).GT.20.) THEN
 
363
               DO 190 J = 1,4
 
364
                  PAR(J) = PARS(J)
 
365
  190          CONTINUE
 
366
               GO TO 3000
 
367
            END IF
 
368
C
 
369
            SQM = 0.
 
370
            DO 200 I = 1,NLIV
 
371
               FG  = PAR(1)*EXP(-4*ALOG(2.)*((I-PAR(2))/PAR(3))**2)
 
372
               SQM = SQM+(RY(I)-FG-PAR(4))**2
 
373
  200       CONTINUE
 
374
C
 
375
            SQM = SQRT(SQM/NLIV)
 
376
            SL  = ABS(SQR-SQM)/SQR
 
377
            SQR = SQM
 
378
            GO TO 3001
 
379
 3000    CONTINUE
 
380
C
 
381
         Y1  = PAR(2)+LY-5
 
382
         FO1 = (FO1+PAR(4))/2.
 
383
         H1  = H1-FO1      
 
384
C
 
385
C *** legge magnitudine v , colori b-v e u-b e nome della
 
386
C     stella prescelta e li scrive sul file f
 
387
 3002    CONTINUE
 
388
         WRITE (STRING,9003) IDNGRP
 
389
         NCH = INDEX(STRING,':')+1               
 
390
         CINPUT = ' '
 
391
         CALL STKPRC(STRING(1:NCH),'INPUTC',1,1,60,IAC,CINPUT,
 
392
     2               KUN,KNUL,ISTAT)
 
393
         ICOM   = MIN(INDEX(CINPUT,',')-1,6)
 
394
         IF (ICOM.LE.0) THEN
 
395
            VV       = 0.0
 
396
            BB       = 0.0
 
397
            UU       = 0.0
 
398
         ELSE
 
399
            ICOM=4
 
400
            CALL GENCNV(CINPUT,2,ICOM,IDUM,RINPUT,DDUM,ISTAT)
 
401
C           CALL USRINP(RINPUT,4,'R',CINPUT)
 
402
C           ICOM = NEL(RINPUT,4)
 
403
            IF (ICOM.GT.0 .AND. ICOM.LT.4) THEN
 
404
               DO 3003 ICOM = ICOM+1, 4
 
405
                  RINPUT(ICOM) = 0.0
 
406
 3003          CONTINUE
 
407
            ENDIF
 
408
            IDNGRP= INT(RINPUT(1))
 
409
            IF (IDNGRP.GT.99999) THEN
 
410
               CALL STTPUT('*** WARNING: Group identification should'//
 
411
     2                     ' be less than 100000; try again ...',ISTAT)
 
412
               GO TO 3002
 
413
 
 
414
            ELSE
 
415
               VV     = RINPUT(2)
 
416
               BB     = RINPUT(3)
 
417
               UU     = RINPUT(4)
 
418
            ENDIF 
 
419
         ENDIF
 
420
C
 
421
C *** prepare the arrays for writting
 
422
         IDNGRP     = IDNGRP
 
423
         IDNCMP(1)  = 101
 
424
         IROW       = IROW + 1
 
425
         B1         = VV+BB
 
426
         U1         = B1+UU
 
427
         LX1        = BEGIN(1)+X1-KDX
 
428
         LY1        = BEGIN(2)+Y1-KDY
 
429
         PARINT(1)  = FLOAT(LX1)
 
430
         PARINT(2)  = FLOAT(LY1)
 
431
         PARINT(3)  = VV
 
432
         PARINT(4)  = B1
 
433
         PARINT(5)  = U1
 
434
         PARINT(6)  = DIMX
 
435
         PARINT(7)  = DIMY
 
436
         PARINT(8)  = 0.0
 
437
         PARINT(9)  = 0.0
 
438
         PARINT(10) = FO1
 
439
         PARINT(11) = 4.0
 
440
         PARINT(12) = 0.0
 
441
         PARINT(13) = 0.0
 
442
         PARINT(14) = 0.0
 
443
         PARINT(15) = NCP
 
444
         PARINT(16) = NHL
 
445
C
 
446
         FITCMP(1)  = H1
 
447
         FITCMP(2)  = FLOAT(KDX)
 
448
         FITCMP(3)  = FLOAT(KDY)
 
449
         FITCMP(4)  = 3.0
 
450
         FITCMP(5)  = 0.0
 
451
         FITCMP(6)  = 0.0
 
452
C           
 
453
         CALL INTWWR(TIDINT,IROW,NCP,NHL)
 
454
         NOBJ       = NOBJ + 1
 
455
         NGRP       = NGRP + 1
 
456
         GO TO 1000
 
457
      ENDIF
 
458
 
 
459
      CALL STSEPI
 
460
      END