Filter.cc 12.3 KB
Newer Older
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1
2
3
4
/*
  This file is part of CDO. CDO is a collection of Operators to
  manipulate and analyse Climate model Data.

Uwe Schulzweida's avatar
Uwe Schulzweida committed
5
  Copyright (C) 2003-2017 Uwe Schulzweida, <uwe.schulzweida AT mpimet.mpg.de>
Uwe Schulzweida's avatar
Uwe Schulzweida committed
6
7
8
9
10
11
12
13
14
15
16
17
18
19
  See COPYING file for copying and redistribution conditions.

  This program is free software; you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation; version 2 of the License.

  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.
*/

/*
   This module contains the following operators:
20
21
22
23
 
      Filter    highpass
      Filter    lowpass
      Filter    bandpass
Uwe Schulzweida's avatar
Uwe Schulzweida committed
24
25
*/

Uwe Schulzweida's avatar
Uwe Schulzweida committed
26
#if defined(HAVE_CONFIG_H)
27
#include "config.h"
Ralf Mueller's avatar
6*8=42    
Ralf Mueller committed
28
29
#endif

30

Ralf Mueller's avatar
Ralf Mueller committed
31
#include <cdi.h>
Uwe Schulzweida's avatar
Uwe Schulzweida committed
32
33
#include "cdo.h"
#include "cdo_int.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
34
#include "statistic.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
35
36
#include "pstream.h"

37
#if defined(HAVE_LIBFFTW3) 
38
39
#include <fftw3.h>
#endif
40

Uwe Schulzweida's avatar
Uwe Schulzweida committed
41

Uwe Schulzweida's avatar
Uwe Schulzweida committed
42
#define  NALLOC_INC  1024
Uwe Schulzweida's avatar
Uwe Schulzweida committed
43
44


Uwe Schulzweida's avatar
Uwe Schulzweida committed
45
/* include from Tinfo.c */
46
void getTimeInc(double jdelta, int vdate0, int vdate1, int *incperiod, int *incunit);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
47

Uwe Schulzweida's avatar
Uwe Schulzweida committed
48
49
static
void create_fmasc(int nts, double fdata, double fmin, double fmax, int *fmasc)
50
51
52
{  
  double dimin = nts*fmin / fdata;
  double dimax = nts*fmax / fdata;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
53

54
55
  int imin = (dimin < 0) ? 0 : (int)floor(dimin);  
  int imax = (ceil(dimax) > nts/2) ? nts/2 : (int) ceil(dimax);
56

57
58
  if ( imin < 0 || imin >= nts ) cdoAbort("Parameter fmin=%g: timestep %d out of bounds (1-%d)!", fmin, imin+1, nts);
  if ( imax < 0 || imax >= nts ) cdoAbort("Parameter fmax=%g: timestep %d out of bounds (1-%d)!", fmax, imax+1, nts);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
59

60
  fmasc[imin] = 1;
61
  for ( int i = imin+1; i <= imax; i++ )  
62
63
    fmasc[i] = fmasc[nts-i] = 1; 
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
64

65
#if defined(HAVE_LIBFFTW3) 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
66
static
67
void filter_fftw(int nts, const int *fmasc, fftw_complex *fft_out, fftw_plan *p_T2S, fftw_plan *p_S2T)
68
{
69
70
  fftw_execute(*p_T2S);

71
  for ( int i = 0; i < nts; i++ )
72
73
    if ( ! fmasc[i] )
      {
74
75
        fft_out[i][0] = 0;
        fft_out[i][1] = 0;
76
      }
77
78
79
80
81
  
  fftw_execute(*p_S2T);
  
  return;
}
82
#endif
83
84
85
86

static
void filter_intrinsic(int nts, const int *fmasc, double *array1, double *array2)
{  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
87
88
  bool lpower2 = ((nts&(nts-1)) == 0);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
89
90
91
92
93
  double *work_r = NULL;
  double *work_i = NULL;

  if ( !lpower2 )
    {
94
95
      work_r = (double*) Malloc(nts*sizeof(double));
      work_i = (double*) Malloc(nts*sizeof(double));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
96
    }
97

Uwe Schulzweida's avatar
Uwe Schulzweida committed
98
99
100
101
  if ( lpower2 )
    fft(array1, array2, nts, 1);
  else
    ft_r(array1, array2, nts, 1, work_r, work_i);
102

103
  for ( int i = 0; i < nts; i++ )
104
105
106
    if ( ! fmasc[i] )
      array1[i] = array2[i] = 0;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
107
108
109
110
111
  if ( lpower2 )
    fft(array1, array2, nts, -1);
  else
    ft_r(array1, array2, nts, -1, work_r, work_i);

112
113
  if ( work_r ) Free(work_r);
  if ( work_i ) Free(work_i);
114
  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
115
  return;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
116
117
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
118

Uwe Schulzweida's avatar
Uwe Schulzweida committed
119
120
void *Filter(void *argument)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
121
  enum {BANDPASS, HIGHPASS, LOWPASS};
122
  const char *tunits[] = {"second", "minute", "hour", "day", "month", "year"};
123
  int iunits[] = {31536000, 525600, 8760, 365, 12, 1};
Uwe Schulzweida's avatar
Uwe Schulzweida committed
124
  int nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
125
  int varID, levelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
126
127
  int nalloc = 0;
  int nmiss;
128
  int incperiod0, incunit0, incunit;
129
  int year0, month0, day0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
130
131
  bool use_fftw = false;
  double fmin = 0, fmax = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
132
  double fdata = 0;
133
  field_type ***vars = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
134
  dtlist_type *dtlist = dtlist_new();
135
136
  typedef struct
  {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
137
138
139
    double *array1;
    double *array2;
#if defined(HAVE_LIBFFTW3) 
140
141
142
143
    fftw_complex *in_fft;
    fftw_complex *out_fft;
    fftw_plan p_T2S;
    fftw_plan p_S2T;
144
#endif
145
146
  } memory_t;
  memory_t *ompmem = NULL;
147
  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
148
149
  cdoInitialize(argument);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
150
151
  cdoOperatorAdd("bandpass",  BANDPASS,  0, NULL);
  cdoOperatorAdd("highpass",  HIGHPASS,  0, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
152
  cdoOperatorAdd("lowpass",   LOWPASS,   0, NULL);
153

Uwe Schulzweida's avatar
Uwe Schulzweida committed
154
155
  int operatorID = cdoOperatorID();
  int operfunc   = cdoOperatorF1(operatorID);
156
157
158
159

  if ( CDO_Use_FFTW )
    {
#if defined(HAVE_LIBFFTW3) 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
160
      if ( cdoVerbose ) cdoPrint("Using fftw3 lib");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
161
      use_fftw = true;
162
#else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
163
      if ( cdoVerbose ) cdoPrint("LIBFFTW3 support not compiled in!");
164
165
#endif
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
166
      
Uwe Schulzweida's avatar
Uwe Schulzweida committed
167
  if ( cdoVerbose && !use_fftw ) cdoPrint("Using intrinsic FFT function!");
168
  
169
  int streamID1 = pstreamOpenRead(cdoStreamName(0));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
170

171
  int vlistID1 = pstreamInqVlist(streamID1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
172
  int vlistID2 = vlistDuplicate(vlistID1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
173

Uwe Schulzweida's avatar
Uwe Schulzweida committed
174
175
  int taxisID1 = vlistInqTaxis(vlistID1);
  int taxisID2 = taxisDuplicate(taxisID1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
176
177
  vlistDefTaxis(vlistID2, taxisID2);

178
  int calendar = taxisInqCalendar(taxisID1);  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
179
 
180
  int nvars = vlistNvars(vlistID1);
181
  
182
  int tsID = 0;    
183
  while ( (nrecs = pstreamInqTimestep(streamID1, tsID)) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
184
185
    {
      if ( tsID >= nalloc )
186
187
        {
          nalloc += NALLOC_INC;
188
          vars   = (field_type ***) Realloc(vars, nalloc*sizeof(field_type **));
189
        }
Cedrick Ansorge's avatar
Cedrick Ansorge committed
190
                       
Uwe Schulzweida's avatar
Uwe Schulzweida committed
191
      dtlist_taxisInqTimestep(dtlist, taxisID1, tsID);
192
   
193
194
      vars[tsID] = field_malloc(vlistID1, FIELD_NONE);
           
195
      for ( int recID = 0; recID < nrecs; recID++ )
196
        {
197
          pstreamInqRecord(streamID1, &varID, &levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
198
199
          int gridID   = vlistInqVarGrid(vlistID1, varID);
          int gridsize = gridInqSize(gridID);
200
          vars[tsID][varID][levelID].ptr = (double*) Malloc(gridsize*sizeof(double));
201
          pstreamReadRecord(streamID1, vars[tsID][varID][levelID].ptr, &nmiss);
202
          vars[tsID][varID][levelID].nmiss = nmiss;
203
          if ( nmiss ) cdoAbort("Missing value support for operators in module Filter not added yet!");
204
205
206
        }

      /* get and check time increment */                   
207
      if ( tsID > 0 )
Cedrick Ansorge's avatar
Cedrick Ansorge committed
208
        {    
209
210
          int incperiod = 0;
          int year, month, day;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
211
212
213
214
	  int vdate0 = dtlist_get_vdate(dtlist, tsID-1);
	  int vdate  = dtlist_get_vdate(dtlist, tsID);
	  int vtime0 = dtlist_get_vtime(dtlist, tsID-1);
	  int vtime  = dtlist_get_vtime(dtlist, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
215

Uwe Schulzweida's avatar
Uwe Schulzweida committed
216
217
          cdiDecodeDate(vdate0, &year0, &month0, &day0);               
          cdiDecodeDate(vdate,  &year,  &month,  &day);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
218

Uwe Schulzweida's avatar
Uwe Schulzweida committed
219
220
221
          juldate_t juldate0 = juldate_encode(calendar, vdate0, vtime0);        
          juldate_t juldate  = juldate_encode(calendar, vdate,  vtime);         
          double jdelta   = juldate_to_seconds(juldate_sub(juldate, juldate0));
222
223
          
          if ( tsID == 1 ) 
Cedrick Ansorge's avatar
Cedrick Ansorge committed
224
            {           
Uwe Schulzweida's avatar
Uwe Schulzweida committed
225
              getTimeInc(jdelta, vdate0, vdate, &incperiod0, &incunit0);
226
              incperiod = incperiod0; 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
227
              if ( incperiod == 0 ) cdoAbort("Time step must be different from zero!");
228
              incunit = incunit0;
229
              if ( cdoVerbose ) cdoPrint("Time step %i %s", incperiod, tunits[incunit]);
230
231
232
              fdata = 1.*iunits[incunit]/incperiod;
            }
          else 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
233
            getTimeInc(jdelta, vdate0, vdate, &incperiod, &incunit);        
Uwe Schulzweida's avatar
Uwe Schulzweida committed
234

Uwe Schulzweida's avatar
Uwe Schulzweida committed
235
236
          if ( calendar != CALENDAR_360DAYS && calendar != CALENDAR_365DAYS && calendar != CALENDAR_366DAYS &&
	       incunit0 < 4 && month == 2 && day == 29 && 
237
238
               ( day0 != day || month0 != month || year0 != year ) )
            {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
239
              cdoWarning("Filtering of multi-year times series doesn't works properly with a standard calendar.");
240
241
              cdoWarning("  Please delete the day %i-02-29 (cdo del29feb)", year);
            }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
242

243
          if ( ! ( incperiod == incperiod0 && incunit == incunit0 ) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
244
            cdoWarning("Time increment in step %i (%d%s) differs from step 1 (%d%s)!",
245
                       tsID, incperiod, tunits[incunit], incperiod0, tunits[incunit0]);        
246
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
247
248
      tsID++;
    }
249
  
250
  int nts = tsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
251
  if ( nts <= 1 ) cdoAbort("Number of time steps <= 1!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
252
253
254

  if ( use_fftw )
    {
255
#if defined(HAVE_LIBFFTW3) 
256
      ompmem = (memory_t*) Malloc(ompNumThreads*sizeof(memory_t));
257
      for ( int i = 0; i < ompNumThreads; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
258
	{
259
260
	  ompmem[i].in_fft  = (fftw_complex*) Malloc(nts*sizeof(fftw_complex));
	  ompmem[i].out_fft = (fftw_complex*) Malloc(nts*sizeof(fftw_complex));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
261
262
263
264
265
266
	  ompmem[i].p_T2S = fftw_plan_dft_1d(nts, ompmem[i].in_fft, ompmem[i].out_fft,  1, FFTW_ESTIMATE);
	  ompmem[i].p_S2T = fftw_plan_dft_1d(nts, ompmem[i].out_fft, ompmem[i].in_fft, -1, FFTW_ESTIMATE);
	}
#endif
    }
  else
267
    {
268
      ompmem = (memory_t*) Malloc(ompNumThreads*sizeof(memory_t));
269
      for ( int i = 0; i < ompNumThreads; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
270
	{
271
272
	  ompmem[i].array1 = (double*) Malloc(nts*sizeof(double));
	  ompmem[i].array2 = (double*) Malloc(nts*sizeof(double));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
273
	}
274
    }
275

Uwe Schulzweida's avatar
Uwe Schulzweida committed
276
277
278
279
  switch(operfunc)
    {
    case BANDPASS: 
      {
280
281
        operatorInputArg("lower and upper bound of frequency band");
        operatorCheckArgc(2);
282
283
        fmin = parameter2double(operatorArgv()[0]);
        fmax = parameter2double(operatorArgv()[1]);
284
        break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
285
286
287
      }
    case HIGHPASS:
      {              
288
289
        operatorInputArg("lower bound of frequency pass");
        operatorCheckArgc(1);
290
        fmin = parameter2double(operatorArgv()[0]);
291
292
        fmax = fdata;
        break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
293
294
      }
    case LOWPASS: 
295
296
297
298
      { 
        operatorInputArg("upper bound of frequency pass");
        operatorCheckArgc(1);
        fmin = 0;
299
        fmax = parameter2double(operatorArgv()[0]);
300
        break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
301
302
      }
    }
303
304

  if ( cdoVerbose ) cdoPrint("fmin=%g  fmax=%g", fmin, fmax);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
305
  
306
  int *fmasc = (int*) Calloc(nts, sizeof(int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
307
  create_fmasc(nts, fdata, fmin, fmax, fmasc);
308

Uwe Schulzweida's avatar
Uwe Schulzweida committed
309
  for ( int varID = 0; varID < nvars; varID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
310
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
311
312
313
      int gridID   = vlistInqVarGrid(vlistID1, varID);
      int gridsize = gridInqSize(gridID);
      int nlevel   = zaxisInqSize(vlistInqVarZaxis(vlistID1, varID));
314
      
Uwe Schulzweida's avatar
Uwe Schulzweida committed
315
      for ( int levelID = 0; levelID < nlevel; levelID++ )
316
        {
317
318
          if ( use_fftw )
            {
319
#if defined(HAVE_LIBFFTW3) 
320
#if defined(_OPENMP)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
321
#pragma omp parallel for default(shared)
322
#endif
323
              for ( int i = 0; i < gridsize; i++ )
324
                {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
325
            	  int ompthID = cdo_omp_get_thread_num();
326

Uwe Schulzweida's avatar
Uwe Schulzweida committed
327
                  for ( int tsID = 0; tsID < nts; tsID++ )                              
328
                    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
329
330
                      ompmem[ompthID].in_fft[tsID][0] = vars[tsID][varID][levelID].ptr[i];
                      ompmem[ompthID].in_fft[tsID][1] = 0;
331
332
                    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
333
                  filter_fftw(nts, fmasc, ompmem[ompthID].out_fft, &ompmem[ompthID].p_T2S, &ompmem[ompthID].p_S2T);
334
                  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
335
336
                  for ( int tsID = 0; tsID < nts; tsID++ )
                    vars[tsID][varID][levelID].ptr[i] = ompmem[ompthID].in_fft[tsID][0] / nts;  
337
338
339
340
341
342
                }
#endif
            }
          else
            {
#if defined(_OPENMP)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
343
#pragma omp parallel for default(shared)
344
#endif
345
              for ( int i = 0; i < gridsize; i++ )  
346
                {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
347
            	  int ompthID = cdo_omp_get_thread_num();
348

Uwe Schulzweida's avatar
Uwe Schulzweida committed
349
                  for ( int tsID = 0; tsID < nts; tsID++ )
350
351
352
353
354
355
                    ompmem[ompthID].array1[tsID] = vars[tsID][varID][levelID].ptr[i];

                  memset(ompmem[ompthID].array2, 0, nts*sizeof(double));

                  filter_intrinsic(nts, fmasc, ompmem[ompthID].array1, ompmem[ompthID].array2);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
356
                  for ( int tsID = 0; tsID < nts; tsID++ )
357
358
359
360
                    vars[tsID][varID][levelID].ptr[i] = ompmem[ompthID].array1[tsID];
                }
            }
        }
361
    }
362

Uwe Schulzweida's avatar
Uwe Schulzweida committed
363
364
365
  if ( use_fftw )
    {
#if defined(HAVE_LIBFFTW3) 
366
      for ( int i = 0; i < ompNumThreads; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
367
	{
368
369
	  Free(ompmem[i].in_fft);
	  Free(ompmem[i].out_fft);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
370
	}
371
      Free(ompmem);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
372
373
374
#endif
    }
  else
375
    {
376
      for ( int i = 0; i < ompNumThreads; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
377
	{
378
379
	  Free(ompmem[i].array1);
	  Free(ompmem[i].array2);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
380
	}
381
      Free(ompmem);
382
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
383

384
  int streamID2 = pstreamOpenWrite(cdoStreamName(1), cdoFiletype());
Uwe Schulzweida's avatar
Uwe Schulzweida committed
385
  
386
  pstreamDefVlist(streamID2, vlistID2);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
387
 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
388
  for ( int tsID = 0; tsID < nts; tsID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
389
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
390
      dtlist_taxisDefTimestep(dtlist, taxisID2, tsID);
391
      pstreamDefTimestep(streamID2, tsID);
392
    
Uwe Schulzweida's avatar
Uwe Schulzweida committed
393
      for ( int varID = 0; varID < nvars; varID++ )
394
        {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
395
396
          int nlevel = zaxisInqSize(vlistInqVarZaxis(vlistID1, varID));
          for ( int levelID = 0; levelID < nlevel; levelID++ )
397
398
399
            {
              if ( vars[tsID][varID][levelID].ptr )
                {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
400
                  int nmiss = vars[tsID][varID][levelID].nmiss;
401
402
                  pstreamDefRecord(streamID2, varID, levelID);
                  pstreamWriteRecord(streamID2, vars[tsID][varID][levelID].ptr, nmiss);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
403

404
                  Free(vars[tsID][varID][levelID].ptr);
405
                  vars[tsID][varID][levelID].ptr = NULL;
406
407
408
                }
            }
        }
409
410

      field_free(vars[tsID], vlistID1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
411
412
    }

413
  if ( vars ) Free(vars);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
414
415

  dtlist_delete(dtlist);
416

417
418
  pstreamClose(streamID2);
  pstreamClose(streamID1);
419

Uwe Schulzweida's avatar
Uwe Schulzweida committed
420
  cdoFinish();
421
  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
422
  return 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
423
}