← Back to branch summary

~ubuntu-branches/ubuntu/feisty/comedilib/feisty

« back to all changes in this revision

Viewing changes to comedi_calibrate/comedi_calibrate.c

  • Committer: Bazaar Package Importer
  • Author(s): David Schleef
  • Date: 2003-09-23 18:11:12 UTC
  • Revision ID: james.westby@ubuntu.com-20030923181112-sat05jyh702rb1at
Tags: upstream-0.7.21
Import upstream version 0.7.21

Show diffs side-by-side

added added

removed removed

Lines of Context:
comedi_calibrate/comedi_calibrate.c
 
1
/*
 
2
   A little auto-calibration utility, for boards
 
3
   that support it.
 
4
 
 
5
   copyright (C) 1999,2000,2001,2002 by David Schleef
 
6
   copyright (C) 2003 by Frank Mori Hess
 
7
 
 
8
 */
 
9
 
 
10
 /***************************************************************************
 
11
 *                                                                         *
 
12
 *   This program is free software; you can redistribute it and/or modify  *
 
13
 *   it under the terms of the GNU Lesser General Public License as        *
 
14
 *   published by                                                          *
 
15
 *   the Free Software Foundation; either version 2.1 of the License, or   *
 
16
 *   (at your option) any later version.                                   *
 
17
 *                                                                         *
 
18
 ***************************************************************************/
 
19
 
 
20
#define _GNU_SOURCE
 
21
 
 
22
#include <stdio.h>
 
23
#include "comedilib.h"
 
24
#include <fcntl.h>
 
25
#include <unistd.h>
 
26
#include <errno.h>
 
27
#include <getopt.h>
 
28
#include <ctype.h>
 
29
#include <math.h>
 
30
#include <stdlib.h>
 
31
#include <string.h>
 
32
#include <assert.h>
 
33
 
 
34
#include "calib.h"
 
35
 
 
36
/* global variables */
 
37
int verbose = 0;
 
38
 
 
39
struct board_struct{
 
40
        char *name;
 
41
        char *id;
 
42
        int (*init_setup)( calibration_setup_t *setup, const char *device_name );
 
43
};
 
44
 
 
45
struct board_struct drivers[] = {
 
46
        { "ni_pcimio",  ni_id,  ni_setup },
 
47
        { "ni_atmio",   ni_id,  ni_setup },
 
48
        { "ni_mio_cs",  ni_id,  ni_setup },
 
49
        { "cb_pcidas",  cb_id,  cb_setup },
 
50
        { "cb_pcidas64",        cb64_id,        cb64_setup },
 
51
        { "ni_labpc", ni_labpc_id,      ni_labpc_setup },
 
52
};
 
53
#define n_drivers (sizeof(drivers)/sizeof(drivers[0]))
 
54
 
 
55
void help(void)
 
56
{
 
57
        printf("comedi_calibrate [options] - autocalibrates a Comedi device\n");
 
58
        printf("  --verbose, -v  \n");
 
59
        printf("  --quiet, -q  \n");
 
60
        printf("  --help, -h  \n");
 
61
        printf("  --file, -f [/dev/comediN] \n");
 
62
        printf("  --save-file, -S [filepath] \n");
 
63
        printf("  --driver-name [driver]  \n");
 
64
        printf("  --device-name [device]  \n");
 
65
        printf("  --[no-]reset  \n");
 
66
        printf("  --[no-]calibrate  \n");
 
67
        printf("  --[no-]dump  \n");
 
68
        printf("  --[no-]results  \n");
 
69
        printf("  --[no-]output  \n");
 
70
}
 
71
 
 
72
typedef struct
 
73
{
 
74
        int verbose;
 
75
        char *file_path;
 
76
        char *save_file_path;
 
77
        char *driver_name;
 
78
        char *device_name;
 
79
        int do_reset;
 
80
        int do_dump;
 
81
        int do_calibrate;
 
82
        int do_results;
 
83
        int do_output;
 
84
        unsigned int subdevice;
 
85
        unsigned int channel;
 
86
        unsigned int range;
 
87
        unsigned int aref;
 
88
} parsed_options_t;
 
89
 
 
90
static void parse_options( int argc, char *argv[], parsed_options_t *settings )
 
91
{
 
92
        int c, index;
 
93
 
 
94
        struct option options[] = {
 
95
                { "verbose", 0, 0, 'v' },
 
96
                { "quiet", 0, 0, 'q' },
 
97
                { "file", 1, 0, 'f' },
 
98
                { "save-file", 1, 0, 'S' },
 
99
                { "help", 0, 0, 'h' },
 
100
                { "driver-name", 1, 0, 0x1000 },
 
101
                { "device-name", 1, 0, 0x1001 },
 
102
                { "reset", 0, &settings->do_reset, 1 },
 
103
                { "no-reset", 0, &settings->do_reset, 0 },
 
104
                { "calibrate", 0, &settings->do_calibrate, 1 },
 
105
                { "no-calibrate", 0, &settings->do_calibrate, 0 },
 
106
                { "dump", 0, &settings->do_dump, 1 },
 
107
                { "no-dump", 0, &settings->do_dump, 0 },
 
108
                { "results", 0, &settings->do_results, 1 },
 
109
                { "no-results", 0, &settings->do_results, 0 },
 
110
                { "output", 0, &settings->do_output, 1 },
 
111
                { "no-output", 0, &settings->do_output, 0 },
 
112
                { "subdevice", 1, 0, 's' },
 
113
                { "channel", 1, 0, 'c' },
 
114
                { "range", 1, 0, 'r' },
 
115
                { "aref", 1, 0, 'a' },
 
116
                { 0 },
 
117
        };
 
118
 
 
119
        while (1) {
 
120
                c = getopt_long(argc, argv, "f:S:vqs:c:r:a:", options, &index);
 
121
                if (c == -1)break;
 
122
                switch (c) {
 
123
                case 0:
 
124
                        continue;
 
125
                case 'h':
 
126
                        help();
 
127
                        exit(0);
 
128
                        break;
 
129
                case 'f':
 
130
                        settings->file_path = optarg;
 
131
                        break;
 
132
                case 'S':
 
133
                        settings->save_file_path = optarg;
 
134
                        break;
 
135
                case 'v':
 
136
                        settings->verbose++;
 
137
                        break;
 
138
                case 'q':
 
139
                        settings->verbose--;
 
140
                        break;
 
141
                case 0x1000:
 
142
                        settings->driver_name = optarg;
 
143
                        break;
 
144
                case 0x1001:
 
145
                        settings->device_name = optarg;
 
146
                        break;
 
147
                case 's':
 
148
                        settings->subdevice = strtoul( optarg, NULL, 0 );
 
149
                        break;
 
150
                case 'c':
 
151
                        settings->channel = strtoul( optarg, NULL, 0 );
 
152
                        break;
 
153
                case 'r':
 
154
                        settings->range = strtoul( optarg, NULL, 0 );
 
155
                        break;
 
156
                case 'a':
 
157
                        settings->aref = strtoul( optarg, NULL, 0 );
 
158
                        break;
 
159
                default:
 
160
                        help();
 
161
                        exit(1);
 
162
                }
 
163
        }
 
164
}
 
165
 
 
166
int main(int argc, char *argv[])
 
167
{
 
168
        int i;
 
169
        struct board_struct *this_board;
 
170
        int device_status = STATUS_UNKNOWN;
 
171
        calibration_setup_t setup;
 
172
        int retval;
 
173
        parsed_options_t options;
 
174
 
 
175
        memset( &setup, 0, sizeof( setup ) );
 
176
        setup.sv_settling_time_ns = 99999;
 
177
        setup.sv_order = 10;
 
178
        
 
179
        memset( &options, 0, sizeof( options ) );
 
180
        options.do_dump = 0;
 
181
        options.do_reset = 0;
 
182
        options.do_calibrate = -1;
 
183
        options.do_results = 0;
 
184
        options.do_output = 1;
 
185
        options.file_path = "/dev/comedi0";
 
186
        parse_options( argc, argv, &options );
 
187
        verbose = options.verbose;
 
188
        
 
189
        setup.dev = comedi_open( options.file_path );
 
190
        if( setup.dev == NULL ) {
 
191
                fprintf( stderr, "comedi_open() failed, with device file name: %s\n",
 
192
                        options.file_path );
 
193
                comedi_perror("comedi_open");
 
194
                exit(0);
 
195
        }
 
196
 
 
197
        if( options.save_file_path == NULL )
 
198
                options.save_file_path = comedi_get_default_calibration_path( setup.dev );
 
199
        if(!options.driver_name)
 
200
                options.driver_name=comedi_get_driver_name( setup.dev );
 
201
        if(!options.device_name)
 
202
                options.device_name=comedi_get_board_name( setup.dev );
 
203
 
 
204
        setup.ad_subdev=comedi_find_subdevice_by_type( setup.dev,COMEDI_SUBD_AI,0);
 
205
        setup.da_subdev=comedi_find_subdevice_by_type( setup.dev,COMEDI_SUBD_AO,0);
 
206
        setup.caldac_subdev=comedi_find_subdevice_by_type( setup.dev,COMEDI_SUBD_CALIB,0);
 
207
        setup.eeprom_subdev=comedi_find_subdevice_by_type( setup.dev,COMEDI_SUBD_MEMORY,0);
 
208
 
 
209
        for(i=0;i<n_drivers;i++){
 
210
                if(!strcmp(drivers[i].name,options.driver_name)){
 
211
                        this_board = drivers+i;
 
212
                        goto ok;
 
213
                }
 
214
        }
 
215
        fprintf(stderr, "Driver %s unknown\n", options.driver_name);
 
216
        return 1;
 
217
 
 
218
ok:
 
219
 
 
220
        retval = this_board->init_setup( &setup, options.device_name );
 
221
        if( retval < 0 ){
 
222
                fprintf(stderr, "init_setup() failed for %s\n", options.device_name );
 
223
                return 1;
 
224
        }
 
225
        device_status = setup.status;
 
226
 
 
227
        if(device_status<STATUS_DONE){
 
228
                printf("Warning: device may not be not fully calibrated due to "
 
229
                        "insufficient information.\n"
 
230
                        "Please file a bug report at https://bugs.comedi.org/ and attach this output.\n"
 
231
                        "This output will also allow comedi_calibrate to execute more\n"
 
232
                        "quickly in the future.\n");
 
233
                if(verbose<1){
 
234
                        verbose=1;
 
235
                        printf("Forcing option: --verbose\n");
 
236
                }
 
237
                if(device_status==STATUS_UNKNOWN){
 
238
                        options.do_reset=1;
 
239
                        options.do_dump=1;
 
240
                        options.do_calibrate=0;
 
241
                        options.do_results=0;
 
242
                        printf("Forcing options: --reset --dump --no-calibrate --no-results\n");
 
243
                }
 
244
                if(device_status==STATUS_SOME){
 
245
                        options.do_reset=1;
 
246
                        options.do_dump=1;
 
247
                        options.do_calibrate=1;
 
248
                        options.do_results=1;
 
249
                        printf("Forcing options: --reset --dump --calibrate --results\n");
 
250
                }
 
251
                if(device_status==STATUS_GUESS){
 
252
                        options.do_reset=1;
 
253
                        options.do_dump=1;
 
254
                        options.do_calibrate=1;
 
255
                        options.do_results=1;
 
256
                        printf("Forcing options: --reset --dump --calibrate --results\n");
 
257
                }
 
258
        }
 
259
        if(verbose>=0){
 
260
                char *s = "$Id: comedi_calibrate.c,v 1.89 2003/08/08 17:38:59 fmhess Exp $";
 
261
 
 
262
                printf("%.*s\n",(int)strlen(s)-2,s+1);
 
263
                printf("Driver name: %s\n", options.driver_name);
 
264
                printf("Device name: %s\n", options.device_name);
 
265
                printf("%.*s\n",(int)strlen(this_board->id)-2,this_board->id+1);
 
266
                printf("Comedi version: %d.%d.%d\n",
 
267
                        (comedi_get_version_code(setup.dev)>>16)&0xff,
 
268
                        (comedi_get_version_code(setup.dev)>>8)&0xff,
 
269
                        (comedi_get_version_code(setup.dev))&0xff);
 
270
        }
 
271
 
 
272
        if( options.do_reset == 0 )
 
273
                setup.old_calibration = comedi_parse_calibration_file( options.save_file_path );
 
274
        else
 
275
                setup.old_calibration = NULL;
 
276
        if( options.do_calibrate < 0 )
 
277
        {
 
278
                if( setup.old_calibration ) options.do_calibrate = 0;
 
279
                else options.do_calibrate = 1;
 
280
        }
 
281
        setup.do_output = options.do_output;
 
282
        setup.cal_save_file_path = options.save_file_path;
 
283
 
 
284
        if(options.do_dump) observe( &setup );
 
285
        if(options.do_calibrate && setup.do_cal)
 
286
        {
 
287
                setup.new_calibration = malloc( sizeof( comedi_calibration_t ) );
 
288
                assert( setup.new_calibration );
 
289
                memset( setup.new_calibration, 0, sizeof( comedi_calibration_t ) );
 
290
                setup.new_calibration->driver_name = strdup( comedi_get_driver_name( setup.dev ) );
 
291
                assert( setup.new_calibration->driver_name != NULL );
 
292
                setup.new_calibration->board_name = strdup( comedi_get_board_name( setup.dev ) );
 
293
                assert( setup.new_calibration->board_name != NULL );
 
294
                retval = setup.do_cal( &setup );
 
295
                if( retval < 0 )
 
296
                {
 
297
                        fprintf( stderr, "calibration function returned error\n" );
 
298
                        return -1;
 
299
                }
 
300
        }
 
301
        if(options.do_results) observe( &setup );
 
302
 
 
303
        if( setup.old_calibration ) comedi_cleanup_calibration( setup.old_calibration );
 
304
        if( setup.new_calibration ) comedi_cleanup_calibration( setup.new_calibration );
 
305
 
 
306
        retval = comedi_apply_calibration( setup.dev, options.subdevice,
 
307
                options.channel, options.range, options.aref, setup.cal_save_file_path );
 
308
        if( retval < 0 )
 
309
        {
 
310
                DPRINT( 0, "Failed to apply " );
 
311
        }else
 
312
        {
 
313
                DPRINT( 0, "Applied " );
 
314
        }
 
315
        DPRINT( 0, "calibration for subdevice %i, channel %i, range %i, aref %i\n",
 
316
                options.subdevice, options.channel, options.range, options.aref );
 
317
        comedi_close(setup.dev);
 
318
 
 
319
        return retval;
 
320
}
 
321
 
 
322
void set_target( calibration_setup_t *setup, int obs,double target)
 
323
{
 
324
        comedi_range *range;
 
325
        lsampl_t maxdata, data;
 
326
 
 
327
        comedi_set_global_oor_behavior( COMEDI_OOR_NUMBER );
 
328
 
 
329
        range = comedi_get_range(setup->dev,
 
330
                setup->observables[obs].preobserve_insn.subdev,
 
331
                CR_CHAN( setup->observables[obs].preobserve_insn.chanspec ),
 
332
                CR_RANGE( setup->observables[obs].preobserve_insn.chanspec ));
 
333
        assert( range );
 
334
        maxdata = comedi_get_maxdata( setup->dev,
 
335
                setup->observables[obs].preobserve_insn.subdev,
 
336
                CR_CHAN(setup->observables[obs].preobserve_insn.chanspec));
 
337
        assert( maxdata > 0 );
 
338
        data = comedi_from_phys(target,range,maxdata);
 
339
 
 
340
        setup->observables[obs].preobserve_data[0] = data;
 
341
        setup->observables[obs].target = comedi_to_phys(data,range,maxdata);
 
342
}
 
343
 
 
344
static void apply_appropriate_cal( calibration_setup_t *setup, comedi_insn insn )
 
345
{
 
346
        int retval = 0;
 
347
 
 
348
        if( setup->new_calibration )
 
349
        {
 
350
                retval = comedi_apply_parsed_calibration( setup->dev, insn.subdev,
 
351
                        CR_CHAN( insn.chanspec ), CR_RANGE( insn.chanspec ),
 
352
                        CR_AREF( insn.chanspec ), setup->new_calibration );
 
353
        }else if( setup->old_calibration )
 
354
        {
 
355
                retval = comedi_apply_parsed_calibration( setup->dev, insn.subdev,
 
356
                        CR_CHAN( insn.chanspec ), CR_RANGE( insn.chanspec ),
 
357
                        CR_AREF( insn.chanspec ), setup->old_calibration );
 
358
        }else
 
359
        {
 
360
                reset_caldacs( setup );
 
361
                return;
 
362
        }
 
363
        if( retval < 0 )
 
364
                DPRINT( 1, "failed to apply ");
 
365
        else
 
366
                DPRINT( 1, "applied ");
 
367
        DPRINT( 1, "calibration for subdev %i, channel %i, range %i, aref %i\n", insn.subdev,
 
368
                CR_CHAN( insn.chanspec ), CR_RANGE( insn.chanspec ),
 
369
                CR_AREF( insn.chanspec ) );
 
370
}
 
371
 
 
372
void observe( calibration_setup_t *setup )
 
373
{
 
374
        int i;
 
375
        observable *obs;
 
376
 
 
377
        for( i = 0; i < setup->n_observables; i++){
 
378
                obs = &setup->observables[ i ];
 
379
                if( obs->observe_insn.n == 0 ) continue;
 
380
                preobserve( setup, i);
 
381
                DPRINT(0,"%s\n", setup->observables[i].name);
 
382
                if( obs->preobserve_insn.n != 0){
 
383
                        apply_appropriate_cal( setup, obs->preobserve_insn );
 
384
                }
 
385
                apply_appropriate_cal( setup, obs->observe_insn );
 
386
                measure_observable( setup, i);
 
387
                if(verbose>=1){
 
388
                        observable_dependence( setup, i);
 
389
                }
 
390
        }
 
391
}
 
392
 
 
393
int preobserve( calibration_setup_t *setup, int obs)
 
394
{
 
395
        int retval;
 
396
        comedi_insn reference_source_config;
 
397
        lsampl_t ref_data[ 2 ];
 
398
 
 
399
        // setup reference source
 
400
        memset( &reference_source_config, 0, sizeof(reference_source_config) );
 
401
        reference_source_config.insn = INSN_CONFIG;
 
402
        reference_source_config.n = 2;
 
403
        reference_source_config.subdev = setup->ad_subdev;
 
404
        reference_source_config.data = ref_data;
 
405
        reference_source_config.data[ 0 ] = INSN_CONFIG_ALT_SOURCE;
 
406
        reference_source_config.data[ 1 ] = setup->observables[obs].reference_source;
 
407
 
 
408
        retval = comedi_do_insn( setup->dev, &reference_source_config );
 
409
        /* ignore errors for now since older ni driver doesn't
 
410
         * support reference config insn */
 
411
        if( retval < 0 )
 
412
                perror("preobserve() ignoring reference config error" );
 
413
        retval = 0;
 
414
 
 
415
        if( setup->observables[obs].preobserve_insn.n != 0){
 
416
                retval = comedi_do_insn( setup->dev, &setup->observables[obs].preobserve_insn);
 
417
        }
 
418
        if( retval < 0 )
 
419
                perror("preobserve()");
 
420
 
 
421
        return retval;
 
422
}
 
423
 
 
424
void measure_observable( calibration_setup_t *setup, int obs)
 
425
{
 
426
        char s[32];
 
427
        int n;
 
428
        new_sv_t sv;
 
429
 
 
430
        my_sv_init(&sv, setup,
 
431
                setup->observables[obs].observe_insn.subdev,
 
432
                setup->observables[obs].observe_insn.chanspec);
 
433
        n = new_sv_measure(setup->dev, &sv);
 
434
 
 
435
        sci_sprint_alt(s,sv.average,sv.error);
 
436
        DPRINT(0,"reading %s, target %g\n",s, setup->observables[obs].target);
 
437
        assert( isnan( setup->observables[obs].target) == 0 );
 
438
}
 
439
 
 
440
void observable_dependence(calibration_setup_t *setup, int obs)
 
441
{
 
442
        int i;
 
443
        linear_fit_t l;
 
444
 
 
445
        for( i = 0; i < setup->n_caldacs; i++){
 
446
                check_gain_chan_x( setup, &l,
 
447
                        setup->observables[obs].observe_insn.chanspec, i);
 
448
        }
 
449
 
 
450
}
 
451
 
 
452
 
 
453
void postgain_cal( calibration_setup_t *setup, int obs1, int obs2, int dac)
 
454
{
 
455
        double offset1,offset2;
 
456
        linear_fit_t l;
 
457
        double slope1,slope2;
 
458
        double a;
 
459
        double gain;
 
460
        comedi_range *range1,*range2;
 
461
 
 
462
        DPRINT(0,"postgain: %s; %s\n", setup->observables[obs1].name,
 
463
                setup->observables[obs2].name);
 
464
        preobserve(setup, obs1);
 
465
        check_gain_chan_x( setup, &l, setup->observables[obs1].observe_insn.chanspec, dac);
 
466
        offset1=linear_fit_func_y(&l, setup->caldacs[dac].current);
 
467
        DPRINT(2,"obs1: [%d] offset %g\n",obs1,offset1);
 
468
        range1 = comedi_get_range(setup->dev, setup->observables[obs1].observe_insn.subdev,
 
469
                CR_CHAN( setup->observables[obs1].observe_insn.chanspec),
 
470
                CR_RANGE( setup->observables[obs1].observe_insn.chanspec));
 
471
        slope1=l.slope;
 
472
 
 
473
        preobserve( setup, obs2);
 
474
        check_gain_chan_x( setup, &l, setup->observables[obs2].observe_insn.chanspec,dac);
 
475
        offset2=linear_fit_func_y(&l, setup->caldacs[dac].current);
 
476
        DPRINT(2,"obs2: [%d] offset %g\n",obs2,offset2);
 
477
        range2 = comedi_get_range(setup->dev, setup->observables[obs2].observe_insn.subdev,
 
478
                CR_CHAN( setup->observables[obs2].observe_insn.chanspec),
 
479
                CR_RANGE( setup->observables[obs2].observe_insn.chanspec));
 
480
        slope2=l.slope;
 
481
 
 
482
        gain = (range1->max-range1->min)/(range2->max-range2->min);
 
483
        DPRINT(4,"range1 %g range2 %g\n", range1->max-range1->min,
 
484
                range2->max-range2->min);
 
485
        DPRINT(3,"gain: %g\n",gain);
 
486
 
 
487
        DPRINT(3,"difference: %g\n",offset2-offset1);
 
488
 
 
489
        a = (offset1-offset2)/(slope1-slope2);
 
490
        a=setup->caldacs[dac].current-a;
 
491
 
 
492
        update_caldac( setup, dac, rint(a) );
 
493
        usleep(caldac_settle_usec);
 
494
 
 
495
        DPRINT(0,"caldac[%d] set to %g (%g)\n",dac,rint(a),a);
 
496
 
 
497
        if(verbose>=2){
 
498
                preobserve( setup, obs1);
 
499
                measure_observable( setup, obs1);
 
500
                preobserve( setup, obs2);
 
501
                measure_observable( setup, obs2);
 
502
        }
 
503
}
 
504
 
 
505
void cal1( calibration_setup_t *setup, int obs, int dac)
 
506
{
 
507
        linear_fit_t l;
 
508
        double a;
 
509
 
 
510
        DPRINT(0,"linear: %s\n", setup->observables[obs].name);
 
511
        preobserve( setup, obs);
 
512
        check_gain_chan_x( setup, &l, setup->observables[obs].observe_insn.chanspec,dac);
 
513
        a=linear_fit_func_x(&l, setup->observables[obs].target);
 
514
 
 
515
        update_caldac( setup, dac, rint(a) );
 
516
        usleep(caldac_settle_usec);
 
517
 
 
518
        DPRINT(0,"caldac[%d] set to %g (%g)\n",dac,rint(a),a);
 
519
        if(verbose>=3){
 
520
                measure_observable( setup, obs);
 
521
        }
 
522
}
 
523
 
 
524
void cal1_fine( calibration_setup_t *setup, int obs, int dac )
 
525
{
 
526
        linear_fit_t l;
 
527
        double a;
 
528
 
 
529
        DPRINT(0,"linear fine: %s\n", setup->observables[obs].name);
 
530
        preobserve( setup, obs);
 
531
        check_gain_chan_fine( setup, &l, setup->observables[obs].observe_insn.chanspec,dac);
 
532
        a=linear_fit_func_x(&l,setup->observables[obs].target);
 
533
 
 
534
        update_caldac( setup, dac, rint(a) );
 
535
        usleep(caldac_settle_usec);
 
536
 
 
537
        DPRINT(0,"caldac[%d] set to %g (%g)\n",dac,rint(a),a);
 
538
        if(verbose>=3){
 
539
                measure_observable( setup, obs);
 
540
        }
 
541
}
 
542
 
 
543
void peg_binary( calibration_setup_t *setup, int obs, int dac, int maximize )
 
544
{
 
545
        int x0, x1, x;
 
546
        double y0, y1;
 
547
        new_sv_t sv;
 
548
        unsigned int chanspec = setup->observables[obs].observe_insn.chanspec;
 
549
        int polarity;
 
550
 
 
551
        DPRINT(0,"binary peg: %s\n", setup->observables[obs].name);
 
552
        preobserve( setup, obs);
 
553
 
 
554
        comedi_set_global_oor_behavior( COMEDI_OOR_NUMBER );
 
555
 
 
556
        my_sv_init(&sv, setup, setup->ad_subdev, chanspec);
 
557
 
 
558
        x0 = setup->caldacs[dac].maxdata;
 
559
        update_caldac( setup, dac, x0 );
 
560
        usleep(caldac_settle_usec);
 
561
        new_sv_measure( setup->dev, &sv);
 
562
        y0 = sv.average;
 
563
 
 
564
        x1 = 0;
 
565
        update_caldac( setup, dac, x1 );
 
566
        usleep(caldac_settle_usec);
 
567
        new_sv_measure( setup->dev, &sv);
 
568
        y1 = sv.average;
 
569
 
 
570
        if( (y0 - y1) > 0.0 ) polarity = 1;
 
571
        else polarity = -1;
 
572
 
 
573
        if( maximize )
 
574
        {
 
575
                if( polarity > 0 ) x = x0;
 
576
                else x = x1;
 
577
        }else
 
578
        {
 
579
                if( polarity > 0 ) x = x1;
 
580
                else x = x0;
 
581
        }
 
582
        update_caldac( setup, dac, x );
 
583
        DPRINT(0,"caldac[%d] set to %d\n",dac,x);
 
584
        if(verbose>=3){
 
585
                measure_observable( setup, obs);
 
586
        }
 
587
}
 
588
 
 
589
void cal_binary( calibration_setup_t *setup, int obs, int dac)
 
590
{
 
591
        int x0, x1, x2, x;
 
592
        unsigned int bit;
 
593
        double y0, y1, y2;
 
594
        new_sv_t sv;
 
595
        double target = setup->observables[obs].target;
 
596
        unsigned int chanspec = setup->observables[obs].observe_insn.chanspec;
 
597
        int polarity;
 
598
 
 
599
        DPRINT(0,"binary: %s\n", setup->observables[obs].name);
 
600
        preobserve( setup, obs);
 
601
 
 
602
        comedi_set_global_oor_behavior( COMEDI_OOR_NUMBER );
 
603
 
 
604
        my_sv_init(&sv, setup, setup->ad_subdev, chanspec);
 
605
 
 
606
        x0 = setup->caldacs[dac].maxdata;
 
607
        update_caldac( setup, dac, x0 );
 
608
        usleep(caldac_settle_usec);
 
609
        new_sv_measure( setup->dev, &sv);
 
610
        y0 = sv.average;
 
611
 
 
612
        x1 = x2 = 0;
 
613
        update_caldac( setup, dac, x1 );
 
614
        usleep(caldac_settle_usec);
 
615
        new_sv_measure( setup->dev, &sv);
 
616
        y1 = y2 = sv.average;
 
617
 
 
618
        if( (y0 - y1) > 0.0 ) polarity = 1;
 
619
        else polarity = -1;
 
620
 
 
621
        bit = 1;
 
622
        while( ( bit << 1 ) < setup->caldacs[dac].maxdata )
 
623
                bit <<= 1;
 
624
        for( ; bit; bit >>= 1 ){
 
625
                x2 = x1 | bit;
 
626
 
 
627
                update_caldac( setup, dac, x2 );
 
628
                usleep(caldac_settle_usec);
 
629
                new_sv_measure( setup->dev, &sv);
 
630
                y2 = sv.average;
 
631
                DPRINT(3,"trying %d, result %g, target %g\n",x2,y2,target);
 
632
 
 
633
                if( (y2 - target) * polarity < 0.0 ){
 
634
                        x1 = x2;
 
635
                        y1 = y2;
 
636
                }
 
637
 
 
638
                if(verbose>=3){
 
639
                        measure_observable( setup, obs);
 
640
                }
 
641
        }
 
642
 
 
643
        // get that least signficant bit right
 
644
        if( fabs( y1 - target ) < fabs( y2 - target ) )
 
645
                x = x1;
 
646
        else x = x2;
 
647
 
 
648
        update_caldac( setup, dac, x );
 
649
        DPRINT(0,"caldac[%d] set to %d\n",dac,x);
 
650
        if( x >= setup->caldacs[dac].maxdata || x <= 0 )
 
651
                DPRINT(0,"WARNING: caldac[%d] pegged!\n", dac );
 
652
        if(verbose>=3){
 
653
                measure_observable( setup, obs);
 
654
        }
 
655
}
 
656
 
 
657
void cal_postgain_binary( calibration_setup_t *setup, int obs1, int obs2, int dac)
 
658
{
 
659
        cal_relative_binary( setup, obs1, obs2, dac );
 
660
}
 
661
 
 
662
void cal_relative_binary( calibration_setup_t *setup, int obs1, int obs2, int dac)
 
663
{
 
664
        int x0, x1, x2, x, polarity;
 
665
        double y0, y1, y2;
 
666
        new_sv_t sv1, sv2;
 
667
        double target = setup->observables[obs1].target - setup->observables[obs2].target;
 
668
        unsigned int chanspec1 = setup->observables[obs1].observe_insn.chanspec;
 
669
        unsigned int chanspec2 = setup->observables[obs2].observe_insn.chanspec;
 
670
        unsigned int bit;
 
671
 
 
672
        DPRINT(0,"relative binary: %s, %s\n", setup->observables[obs1].name,
 
673
                setup->observables[obs2].name);
 
674
 
 
675
        comedi_set_global_oor_behavior( COMEDI_OOR_NUMBER );
 
676
 
 
677
        x0 = setup->caldacs[dac].maxdata;
 
678
        update_caldac( setup, dac, x0 );
 
679
        usleep(caldac_settle_usec);
 
680
        preobserve( setup, obs1);
 
681
        my_sv_init(&sv1, setup, setup->ad_subdev,chanspec1);
 
682
        new_sv_measure( setup->dev, &sv1);
 
683
 
 
684
        preobserve( setup, obs2);
 
685
        my_sv_init(&sv2, setup, setup->ad_subdev,chanspec2);
 
686
        new_sv_measure( setup->dev, &sv2);
 
687
        y0 = sv1.average - sv2.average;
 
688
 
 
689
        x1 = x2 = 0;
 
690
        update_caldac( setup, dac, x1 );
 
691
        usleep(caldac_settle_usec);
 
692
        preobserve( setup, obs1);
 
693
        new_sv_measure( setup->dev, &sv1);
 
694
 
 
695
        preobserve( setup, obs2);
 
696
        new_sv_measure( setup->dev, &sv2);
 
697
        y1 = y2 = sv1.average - sv2.average;
 
698
 
 
699
        if( (y0 - y1) > 0.0 ) polarity = 1;
 
700
        else polarity = -1;
 
701
 
 
702
        bit = 1;
 
703
        while( ( bit << 1 ) < setup->caldacs[dac].maxdata )
 
704
                bit <<= 1;
 
705
        for( ; bit; bit >>= 1 )
 
706
        {
 
707
                x2 = x1 | bit;
 
708
 
 
709
                update_caldac( setup, dac, x2 );
 
710
                usleep(caldac_settle_usec);
 
711
 
 
712
                preobserve( setup, obs1);
 
713
                new_sv_measure( setup->dev, &sv1);
 
714
 
 
715
                preobserve( setup, obs2);
 
716
                new_sv_measure( setup->dev, &sv2);
 
717
                y2 = sv1.average - sv2.average;
 
718
 
 
719
                DPRINT(3,"trying %d, result %g, target %g\n",x2,y2,target);
 
720
 
 
721
                if( (y2 - target) * polarity < 0.0 ){
 
722
                        x1 = x2;
 
723
                        y1 = y2;
 
724
                }
 
725
 
 
726
                if(verbose>=3){
 
727
                        preobserve( setup, obs1);
 
728
                        measure_observable( setup, obs1);
 
729
                        preobserve( setup, obs2);
 
730
                        measure_observable( setup, obs2);
 
731
                }
 
732
        }
 
733
 
 
734
        if( fabs( y1 - target ) < fabs( y2 - target ) )
 
735
                x = x1;
 
736
        else
 
737
                x = x2;
 
738
        update_caldac( setup, dac, x );
 
739
        DPRINT(0,"caldac[%d] set to %d\n",dac,x);
 
740
        if( x >= setup->caldacs[dac].maxdata || x <= 0 )
 
741
                DPRINT(0,"WARNING: caldac[%d] pegged!\n", dac );
 
742
        if(verbose>=3){
 
743
                preobserve( setup, obs1);
 
744
                measure_observable( setup, obs1);
 
745
                preobserve( setup, obs2);
 
746
                measure_observable( setup, obs2);
 
747
        }
 
748
}
 
749
 
 
750
void cal_linearity_binary( calibration_setup_t *setup, int obs1, int obs2, int obs3, int dac)
 
751
{
 
752
        int x0, x1, x2, x, polarity;
 
753
        double y0, y1, y2;
 
754
        new_sv_t sv1, sv2, sv3;
 
755
        double target = ( setup->observables[obs3].target - setup->observables[obs2].target ) /
 
756
                ( setup->observables[obs2].target - setup->observables[obs1].target );
 
757
        unsigned int chanspec1 = setup->observables[obs1].observe_insn.chanspec;
 
758
        unsigned int chanspec2 = setup->observables[obs2].observe_insn.chanspec;
 
759
        unsigned int chanspec3 = setup->observables[obs3].observe_insn.chanspec;
 
760
        unsigned int bit;
 
761
 
 
762
        DPRINT(0,"linearity binary: %s,\n%s,\n%s\n", setup->observables[obs1].name,
 
763
                setup->observables[obs2].name,setup->observables[obs3].name);
 
764
 
 
765
        comedi_set_global_oor_behavior( COMEDI_OOR_NUMBER );
 
766
 
 
767
        x0 = setup->caldacs[dac].maxdata;
 
768
        update_caldac( setup, dac, x0 );
 
769
        usleep(caldac_settle_usec);
 
770
 
 
771
        preobserve( setup, obs1);
 
772
        my_sv_init(&sv1, setup, setup->ad_subdev,chanspec1);
 
773
        new_sv_measure( setup->dev, &sv1);
 
774
 
 
775
        preobserve( setup, obs2);
 
776
        my_sv_init(&sv2, setup, setup->ad_subdev,chanspec2);
 
777
        new_sv_measure( setup->dev, &sv2);
 
778
 
 
779
        preobserve( setup, obs3);
 
780
        my_sv_init(&sv3, setup, setup->ad_subdev,chanspec3);
 
781
        new_sv_measure( setup->dev, &sv3);
 
782
 
 
783
        y0 = ( sv3.average - sv2.average ) / ( sv2.average - sv1.average );
 
784
 
 
785
        x1 = x2 = 0;
 
786
        update_caldac( setup, dac, x1 );
 
787
        usleep(caldac_settle_usec);
 
788
 
 
789
        preobserve( setup, obs1);
 
790
        new_sv_measure( setup->dev, &sv1);
 
791
 
 
792
        preobserve( setup, obs2);
 
793
        new_sv_measure( setup->dev, &sv2);
 
794
 
 
795
        preobserve( setup, obs3);
 
796
        new_sv_measure( setup->dev, &sv3);
 
797
 
 
798
        y1 = y2 = ( sv3.average - sv2.average ) / ( sv2.average - sv1.average );
 
799
 
 
800
        if( (y0 - y1) > 0.0 ) polarity = 1;
 
801
        else polarity = -1;
 
802
 
 
803
        bit = 1;
 
804
        while( ( bit << 1 ) < setup->caldacs[dac].maxdata )
 
805
                bit <<= 1;
 
806
        for( ; bit; bit >>= 1 )
 
807
        {
 
808
                x2 = x1 | bit;
 
809
 
 
810
                update_caldac( setup, dac, x2 );
 
811
                usleep(caldac_settle_usec);
 
812
 
 
813
                preobserve( setup, obs1);
 
814
                new_sv_measure( setup->dev, &sv1);
 
815
 
 
816
                preobserve( setup, obs2);
 
817
                new_sv_measure( setup->dev, &sv2);
 
818
 
 
819
                preobserve( setup, obs3);
 
820
                new_sv_measure( setup->dev, &sv3);
 
821
 
 
822
                y2 = ( sv3.average - sv2.average ) / ( sv2.average - sv1.average );
 
823
 
 
824
                DPRINT(3,"trying %d, result %g, target %g\n",x2,y2,target);
 
825
 
 
826
                if( (y2 - target) * polarity < 0.0 ){
 
827
                        x1 = x2;
 
828
                        y1 = y2;
 
829
                }
 
830
 
 
831
                if(verbose>=3){
 
832
                        preobserve( setup, obs1);
 
833
                        measure_observable( setup, obs1);
 
834
                        preobserve( setup, obs2);
 
835
                        measure_observable( setup, obs2);
 
836
                        preobserve( setup, obs3);
 
837
                        measure_observable( setup, obs3);
 
838
                }
 
839
        }
 
840
 
 
841
        if( fabs( y1 - target ) < fabs( y2 - target ) )
 
842
                x = x1;
 
843
        else
 
844
                x = x2;
 
845
        update_caldac( setup, dac, x );
 
846
        DPRINT(0,"caldac[%d] set to %d\n",dac,x);
 
847
        if( x >= setup->caldacs[dac].maxdata || x <= 0 )
 
848
                DPRINT(0,"WARNING: caldac[%d] pegged!\n", dac );
 
849
        if(verbose>=3){
 
850
                preobserve( setup, obs1);
 
851
                measure_observable( setup, obs1);
 
852
                preobserve( setup, obs2);
 
853
                measure_observable( setup, obs2);
 
854
                preobserve( setup, obs3);
 
855
                measure_observable( setup, obs3);
 
856
        }
 
857
}
 
858
 
 
859
#if 0
 
860
void chan_cal(int adc,int cdac,int range,double target)
 
861
{
 
862
        linear_fit_t l;
 
863
        double offset;
 
864
        double gain;
 
865
        double a;
 
866
        char s[32];
 
867
 
 
868
        check_gain_chan_x(&l,CR_PACK(adc,range,AREF_OTHER),cdac);
 
869
        offset=linear_fit_func_y(&l,caldacs[cdac].current);
 
870
        gain=l.slope;
 
871
 
 
872
        a=caldacs[cdac].current+(target-offset)/gain;
 
873
 
 
874
        update_caldac( setup, cdac, rint(a));
 
875
 
 
876
        read_chan2(s,adc,range);
 
877
        DPRINT(1,"caldac[%d] set to %g, offset=%s\n",cdac,a,s);
 
878
}
 
879
#endif
 
880
 
 
881
 
 
882
#if 0
 
883
void channel_dependence(int adc,int range)
 
884
{
 
885
        int i;
 
886
        double gain;
 
887
 
 
888
        for(i=0;i<n_caldacs;i++){
 
889
                gain=check_gain_chan(adc,range,i);
 
890
        }
 
891
}
 
892
#endif
 
893
 
 
894
#if 0
 
895
void caldac_dependence(int caldac)
 
896
{
 
897
        int i;
 
898
        double gain;
 
899
 
 
900
        for(i=0;i<8;i++){
 
901
                gain=check_gain_chan(i,0,caldac);
 
902
        }
 
903
}
 
904
#endif
 
905
 
 
906
 
 
907
void setup_caldacs( calibration_setup_t *setup, int caldac_subdev )
 
908
{
 
909
        int n_chan, i;
 
910
 
 
911
        if( caldac_subdev < 0 ){
 
912
                printf("no calibration subdevice\n");
 
913
                return;
 
914
        }
 
915
 
 
916
        // XXX check subdevice type is really calibration
 
917
        // XXX check we dont exceed max number of allowable caldacs
 
918
 
 
919
        n_chan = comedi_get_n_channels( setup->dev, caldac_subdev );
 
920
 
 
921
        for(i = 0; i < n_chan; i++){
 
922
                setup->caldacs[ setup->n_caldacs + i ].subdev = caldac_subdev;
 
923
                setup->caldacs[ setup->n_caldacs + i ].chan = i;
 
924
                setup->caldacs[ setup->n_caldacs + i ].maxdata = comedi_get_maxdata( setup->dev, caldac_subdev, i);
 
925
                setup->caldacs[ setup->n_caldacs + i ].current=0;
 
926
        }
 
927
 
 
928
        setup->n_caldacs += n_chan;
 
929
}
 
930
 
 
931
void reset_caldac( calibration_setup_t *setup, int caldac_index )
 
932
{
 
933
        if( caldac_index < 0 ) return;
 
934
        assert( caldac_index < setup->n_caldacs );
 
935
        update_caldac( setup, caldac_index, setup->caldacs[ caldac_index ].maxdata / 2 );
 
936
}
 
937
 
 
938
void reset_caldacs( calibration_setup_t *setup )
 
939
{
 
940
        int i;
 
941
 
 
942
        for( i = 0; i < setup->n_caldacs; i++){
 
943
                reset_caldac( setup, i );
 
944
        }
 
945
}
 
946
 
 
947
void update_caldac( calibration_setup_t *setup, int caldac_index,
 
948
        int value )
 
949
{
 
950
        int ret;
 
951
        caldac_t *dac;
 
952
 
 
953
        if( caldac_index < 0 ) return;
 
954
        if( caldac_index > setup->n_caldacs )
 
955
        {
 
956
                fprintf( stderr, "invalid caldac index\n" );
 
957
                return;
 
958
        }
 
959
        dac = &setup->caldacs[ caldac_index ];
 
960
        dac->current = value;
 
961
        DPRINT(4,"update %d %d %d\n", dac->subdev, dac->chan, dac->current);
 
962
        if( dac->current < 0 ){
 
963
                DPRINT(1,"caldac set out of range (%d<0)\n", dac->current);
 
964
                dac->current = 0;
 
965
        }
 
966
        if( dac->current > dac->maxdata ){
 
967
                DPRINT(1,"caldac set out of range (%d>%d)\n",
 
968
                        dac->current, dac->maxdata);
 
969
                dac->current = dac->maxdata;
 
970
        }
 
971
 
 
972
        ret = comedi_data_write( setup->dev, dac->subdev, dac->chan, 0, 0,
 
973
                dac->current);
 
974
        if(ret < 0) perror("update_caldac()");
 
975
}
 
976
 
 
977
#if 0
 
978
void check_gain(int ad_chan,int range)
 
979
{
 
980
        int i;
 
981
 
 
982
        for(i=0;i<n_caldacs;i++){
 
983
                check_gain_chan(ad_chan,range,i);
 
984
        }
 
985
}
 
986
#endif
 
987
 
 
988
#if 0
 
989
double check_gain_chan(int ad_chan,int range,int cdac)
 
990
{
 
991
        linear_fit_t l;
 
992
 
 
993
        return check_gain_chan_x(&l,CR_PACK(ad_chan,range,AREF_OTHER),cdac);
 
994
}
 
995
#endif
 
996
 
 
997
 
 
998
double check_gain_chan_x( calibration_setup_t *setup, linear_fit_t *l,unsigned int ad_chanspec,int cdac)
 
999
{
 
1000
        int orig,i,n;
 
1001
        int step;
 
1002
        new_sv_t sv;
 
1003
        double sum_err;
 
1004
        int sum_err_count=0;
 
1005
        char str[20];
 
1006
 
 
1007
        n = setup->caldacs[cdac].maxdata+1;
 
1008
        memset(l,0,sizeof(*l));
 
1009
 
 
1010
        step=n/16;
 
1011
        if(step<1)step=1;
 
1012
        l->n=0;
 
1013
 
 
1014
        l->y_data=malloc(n*sizeof(double)/step);
 
1015
        if(l->y_data == NULL)
 
1016
        {
 
1017
                perror( __FUNCTION__ );
 
1018
                exit(1);
 
1019
        }
 
1020
 
 
1021
        orig = setup->caldacs[cdac].current;
 
1022
 
 
1023
        my_sv_init(&sv, setup, setup->ad_subdev,ad_chanspec);
 
1024
 
 
1025
        update_caldac( setup, cdac, 0 );
 
1026
        usleep(caldac_settle_usec);
 
1027
 
 
1028
        new_sv_measure( setup->dev, &sv);
 
1029
 
 
1030
        sum_err=0;
 
1031
        for(i=0;i*step<n;i++){
 
1032
                update_caldac( setup, cdac, i*step );
 
1033
                usleep(caldac_settle_usec);
 
1034
 
 
1035
                new_sv_measure( setup->dev, &sv);
 
1036
 
 
1037
                l->y_data[i]=sv.average;
 
1038
                if(!isnan(sv.average)){
 
1039
                        sum_err+=sv.error;
 
1040
                        sum_err_count++;
 
1041
                }
 
1042
                l->n++;
 
1043
        }
 
1044
 
 
1045
        update_caldac( setup, cdac, orig );
 
1046
 
 
1047
        l->yerr=sum_err/sqrt(sum_err_count);
 
1048
        l->dx=step;
 
1049
        l->x0=0;
 
1050
 
 
1051
        linear_fit_monotonic(l);
 
1052
 
 
1053
        if(verbose>=2 || (verbose>=1 && fabs(l->slope/l->err_slope)>4.0)){
 
1054
                sci_sprint_alt(str,l->slope,l->err_slope);
 
1055
                printf("caldac[%d] gain=%s V/bit S_min=%g dof=%g\n",
 
1056
                        cdac,str,l->S_min,l->dof);
 
1057
                //printf("--> %g\n",fabs(l.slope/l.err_slope));
 
1058
        }
 
1059
 
 
1060
        if(verbose>=3)dump_curve(l);
 
1061
 
 
1062
        free(l->y_data);
 
1063
 
 
1064
        return l->slope;
 
1065
}
 
1066
 
 
1067
 
 
1068
double check_gain_chan_fine( calibration_setup_t *setup, linear_fit_t *l,unsigned int ad_chanspec,int cdac)
 
1069
{
 
1070
        int orig,i,n;
 
1071
        int step;
 
1072
        new_sv_t sv;
 
1073
        double sum_err;
 
1074
        int sum_err_count=0;
 
1075
        char str[20];
 
1076
        int fine_size = 10;
 
1077
 
 
1078
        n=2*fine_size+1;
 
1079
        memset(l,0,sizeof(*l));
 
1080
 
 
1081
        step=1;
 
1082
        l->n=0;
 
1083
 
 
1084
        l->y_data=malloc(n*sizeof(double)/step);
 
1085
        if(l->y_data == NULL)
 
1086
        {
 
1087
                perror( __FUNCTION__ );
 
1088
                exit(1);
 
1089
        }
 
1090
 
 
1091
        orig = setup->caldacs[cdac].current;
 
1092
 
 
1093
        my_sv_init(&sv, setup, setup->ad_subdev,ad_chanspec);
 
1094
 
 
1095
        update_caldac( setup, cdac, 0 );
 
1096
        usleep(caldac_settle_usec);
 
1097
 
 
1098
        new_sv_measure( setup->dev, &sv);
 
1099
 
 
1100
        sum_err=0;
 
1101
        for(i=0;i<n;i++){
 
1102
                update_caldac( setup, cdac, i+orig-fine_size );
 
1103
                usleep(caldac_settle_usec);
 
1104
 
 
1105
                new_sv_measure( setup->dev, &sv);
 
1106
 
 
1107
                l->y_data[i]=sv.average;
 
1108
                if(!isnan(sv.average)){
 
1109
                        sum_err+=sv.error;
 
1110
                        sum_err_count++;
 
1111
                }
 
1112
                l->n++;
 
1113
        }
 
1114
 
 
1115
        update_caldac( setup, cdac, orig );
 
1116
 
 
1117
        l->yerr=sum_err/sqrt(sum_err_count);
 
1118
        l->dx=1;
 
1119
        l->x0=orig-fine_size;
 
1120
 
 
1121
        linear_fit_monotonic(l);
 
1122
 
 
1123
        if(verbose>=2 || (verbose>=1 && fabs(l->slope/l->err_slope)>4.0)){
 
1124
                sci_sprint_alt(str,l->slope,l->err_slope);
 
1125
                printf("caldac[%d] gain=%s V/bit S_min=%g dof=%g\n",
 
1126
                        cdac,str,l->S_min,l->dof);
 
1127
                //printf("--> %g\n",fabs(l.slope/l.err_slope));
 
1128
        }
 
1129
 
 
1130
        if(verbose>=3)dump_curve(l);
 
1131
 
 
1132
        free(l->y_data);
 
1133
 
 
1134
        return l->slope;
 
1135
}
 
1136
 
 
1137
 
 
1138
 
 
1139
/* helpers */
 
1140
 
 
1141
int is_unipolar( comedi_t *dev, unsigned int subdevice,
 
1142
        unsigned int channel, unsigned int range )
 
1143
{
 
1144
        comedi_range *range_ptr;
 
1145
 
 
1146
        range_ptr = comedi_get_range( dev, subdevice, channel, range );
 
1147
        assert( range_ptr != NULL );
 
1148
        /* This method is better than a direct test, which might fail */
 
1149
        if( fabs( range_ptr->min ) < fabs( range_ptr->max * 0.001 ) )
 
1150
                return 1;
 
1151
        else
 
1152
                return 0;
 
1153
}
 
1154
 
 
1155
int is_bipolar( comedi_t *dev, unsigned int subdevice,
 
1156
        unsigned int channel, unsigned int range )
 
1157
{
 
1158
        comedi_range *range_ptr;
 
1159
 
 
1160
        range_ptr = comedi_get_range( dev, subdevice, channel, range );
 
1161
        assert( range_ptr != NULL );
 
1162
        /* This method is better than a direct test, which might fail */
 
1163
        if( fabs( range_ptr->max + range_ptr->min ) < fabs( range_ptr->max * 0.001 ) )
 
1164
                return 1;
 
1165
        else
 
1166
                return 0;
 
1167
}
 
1168
 
 
1169
int get_bipolar_lowgain(comedi_t *dev,int subdev)
 
1170
{
 
1171
        int ret = -1;
 
1172
        int i;
 
1173
        int n_ranges = comedi_get_n_ranges(dev,subdev,0);
 
1174
        double max = 0;
 
1175
        comedi_range *range;
 
1176
 
 
1177
        for(i=0;i<n_ranges;i++){
 
1178
                range = comedi_get_range(dev,subdev,0,i);
 
1179
                if( is_bipolar( dev, subdev, 0, i ) == 0 ) continue;
 
1180
                if(range->max>max){
 
1181
                        ret = i;
 
1182
                        max=range->max;
 
1183
                }
 
1184
        }
 
1185
 
 
1186
        return ret;
 
1187
}
 
1188
 
 
1189
int get_bipolar_highgain(comedi_t *dev,int subdev)
 
1190
{
 
1191
        int ret = -1;
 
1192
        int i;
 
1193
        int n_ranges = comedi_get_n_ranges(dev,subdev,0);
 
1194
        double min = HUGE_VAL;
 
1195
        comedi_range *range;
 
1196
 
 
1197
        for(i=0;i<n_ranges;i++){
 
1198
                range = comedi_get_range(dev,subdev,0,i);
 
1199
                if( is_bipolar( dev, subdev, 0, i ) == 0 ) continue;
 
1200
                if(range->max<min){
 
1201
                        ret = i;
 
1202
                        min=range->max;
 
1203
                }
 
1204
        }
 
1205
 
 
1206
        return ret;
 
1207
}
 
1208
 
 
1209
int get_unipolar_lowgain(comedi_t *dev,int subdev)
 
1210
{
 
1211
        int ret = -1;
 
1212
        int i;
 
1213
        int n_ranges = comedi_get_n_ranges(dev,subdev,0);
 
1214
        double max = 0.0;
 
1215
        comedi_range *range;
 
1216
 
 
1217
        for(i=0;i<n_ranges;i++){
 
1218
                range = comedi_get_range(dev,subdev,0,i);
 
1219
                if( is_unipolar( dev, subdev, 0, i ) == 0 ) continue;
 
1220
                if(range->max>max){
 
1221
                        ret = i;
 
1222
                        max=range->max;
 
1223
                }
 
1224
        }
 
1225
 
 
1226
        return ret;
 
1227
}
 
1228
 
 
1229
int get_unipolar_highgain(comedi_t *dev,int subdev)
 
1230
{
 
1231
        int ret = -1;
 
1232
        int i;
 
1233
        int n_ranges = comedi_get_n_ranges(dev,subdev,0);
 
1234
        double max = HUGE_VAL;
 
1235
        comedi_range *range;
 
1236
 
 
1237
        for(i=0;i<n_ranges;i++){
 
1238
                range = comedi_get_range(dev,subdev,0,i);
 
1239
                if( is_unipolar( dev, subdev, 0, i ) == 0 ) continue;
 
1240
                if(range->max < max){
 
1241
                        ret = i;
 
1242
                        max=range->max;
 
1243
                }
 
1244
        }
 
1245
 
 
1246
        return ret;
 
1247
}
 
1248
 
 
1249
int read_eeprom( calibration_setup_t *setup, int addr)
 
1250
{
 
1251
        lsampl_t data = 0;
 
1252
        int retval;
 
1253
 
 
1254
        retval = comedi_data_read( setup->dev, setup->eeprom_subdev, addr,0,0,&data);
 
1255
        if( retval < 0 )
 
1256
        {
 
1257
                perror( "read_eeprom()" );
 
1258
                return retval;
 
1259
        }
 
1260
 
 
1261
        return data;
 
1262
}
 
1263
 
 
1264
double read_chan( calibration_setup_t *setup, int adc,int range)
 
1265
{
 
1266
        int n;
 
1267
        new_sv_t sv;
 
1268
        char str[20];
 
1269
 
 
1270
        my_sv_init(&sv, setup,  setup->ad_subdev,CR_PACK(adc,range,AREF_OTHER));
 
1271
 
 
1272
        n=new_sv_measure( setup->dev, &sv);
 
1273
 
 
1274
        sci_sprint_alt(str,sv.average,sv.error);
 
1275
        printf("chan=%d ave=%s\n",adc,str);
 
1276
 
 
1277
        return sv.average;
 
1278
}
 
1279
 
 
1280
int read_chan2( calibration_setup_t *setup, char *s,int adc,int range)
 
1281
{
 
1282
        int n;
 
1283
        new_sv_t sv;
 
1284
 
 
1285
        my_sv_init(&sv, setup, setup->ad_subdev,CR_PACK(adc,range,AREF_OTHER));
 
1286
 
 
1287
        n=new_sv_measure( setup->dev, &sv);
 
1288
 
 
1289
        return sci_sprint_alt(s,sv.average,sv.error);
 
1290
}
 
1291
 
 
1292
#if 0
 
1293
void set_ao(comedi_t *dev,int subdev,int chan,int range,double value)
 
1294
{
 
1295
        comedi_range *r = comedi_get_range(dev,subdev,chan,range);
 
1296
        lsampl_t maxdata = comedi_get_maxdata(dev,subdev,chan);
 
1297
        lsampl_t data;
 
1298
 
 
1299
        data = comedi_from_phys(value,r,maxdata);
 
1300
 
 
1301
        comedi_data_write(dev,subdev,chan,range,AREF_GROUND,data);
 
1302
}
 
1303
#endif
 
1304
 
 
1305
int my_sv_init( new_sv_t *sv, const calibration_setup_t *setup, int subdev,
 
1306
        unsigned int chanspec )
 
1307
{
 
1308
        int retval;
 
1309
        retval = new_sv_init( sv, setup->dev, subdev, chanspec );
 
1310
        sv->settling_time_ns = setup->sv_settling_time_ns;
 
1311
        sv->order = setup->sv_order;
 
1312
        return retval;
 
1313
}
 
1314
 
 
1315
int new_sv_init(new_sv_t *sv,comedi_t *dev,int subdev,unsigned int chanspec)
 
1316
{
 
1317
        memset(sv,0,sizeof(*sv));
 
1318
 
 
1319
        sv->subd=subdev;
 
1320
        //sv->t.flags=TRIG_DITHER;
 
1321
        sv->chanspec = chanspec;
 
1322
 
 
1323
        //sv->chanlist[0]=CR_PACK(chan,range,aref);
 
1324
 
 
1325
        sv->maxdata=comedi_get_maxdata(dev,subdev,CR_CHAN(chanspec));
 
1326
        sv->rng=comedi_get_range(dev,subdev,
 
1327
                CR_CHAN(chanspec), CR_RANGE(chanspec));
 
1328
 
 
1329
        sv->order=10;
 
1330
 
 
1331
        return 0;
 
1332
}
 
1333
 
 
1334
int new_sv_measure( comedi_t *dev, new_sv_t *sv)
 
1335
{
 
1336
        lsampl_t *data;
 
1337
        int n,i,ret;
 
1338
 
 
1339
        double x,s,s2;
 
1340
 
 
1341
        n=1<<sv->order;
 
1342
 
 
1343
        data=malloc(sizeof(lsampl_t)*n);
 
1344
        if(data == NULL)
 
1345
        {
 
1346
                perror( __FUNCTION__ );
 
1347
                exit(1);
 
1348
        }
 
1349
 
 
1350
        ret = comedi_data_read_hint(dev, sv->subd, sv->chanspec, 0, 0);
 
1351
        if(ret<0){
 
1352
                printf("hint barf\n");
 
1353
                goto out;
 
1354
        }
 
1355
        comedi_nanodelay(dev, sv->settling_time_ns);
 
1356
 
 
1357
        ret = comedi_data_read_n(dev, sv->subd, sv->chanspec, 0, 0, data, n);
 
1358
        if(ret<0){
 
1359
                printf("barf\n");
 
1360
                goto out;
 
1361
        }
 
1362
 
 
1363
        s=0.0;
 
1364
        s2=0.0;
 
1365
        for(i = 0; i < n; i++){
 
1366
                x = comedi_to_phys(data[i], sv->rng, sv->maxdata);
 
1367
                s += x;
 
1368
                s2 += x * x;
 
1369
        }
 
1370
        s /= n;
 
1371
        s2 /= n;
 
1372
        sv->average=s;
 
1373
        sv->stddev=sqrt( ( ( n + 1 ) / n ) * ( s2 - s * s ) );
 
1374
        sv->error=sv->stddev / sqrt( n );
 
1375
 
 
1376
        ret=n;
 
1377
 
 
1378
out:
 
1379
        free(data);
 
1380
 
 
1381
        return ret;
 
1382
}
 
1383
 
 
1384
int new_sv_measure_order( comedi_t *dev, new_sv_t *sv,int order)
 
1385
{
 
1386
        lsampl_t *data;
 
1387
        int n,i,ret;
 
1388
        double x,s,s2;
 
1389
 
 
1390
        n=1<<order;
 
1391
 
 
1392
        data=malloc(sizeof(lsampl_t)*n);
 
1393
        if(data == NULL)
 
1394
        {
 
1395
                perror( __FUNCTION__ );
 
1396
                exit(1);
 
1397
        }
 
1398
 
 
1399
        ret = comedi_data_read_n(dev, sv->subd, sv->chanspec, 0, 0, data, n);
 
1400
        if(ret<0){
 
1401
                printf("barf order\n");
 
1402
                goto out;
 
1403
        }
 
1404
 
 
1405
        s=0;
 
1406
        s2=0;
 
1407
        for(i = 0; i < n; i++){
 
1408
                x = comedi_to_phys(data[i], sv->rng, sv->maxdata);
 
1409
                s += x;
 
1410
                s2 += x * x;
 
1411
        }
 
1412
        s /= n;
 
1413
        s2 /= n;
 
1414
        sv->average = s;
 
1415
        sv->stddev=sqrt( ( ( n + 1 ) / n ) * ( s2 - s * s ) );
 
1416
        sv->error=sv->stddev / sqrt( n );
 
1417
 
 
1418
        ret=n;
 
1419
 
 
1420
out:
 
1421
        free(data);
 
1422
 
 
1423
        return ret;
 
1424
}
 
1425
 
 
1426
 
 
1427
 
 
1428
/* linear fitting */
 
1429
 
 
1430
int calculate_residuals(linear_fit_t *l);
 
1431
 
 
1432
int linear_fit_monotonic(linear_fit_t *l)
 
1433
{
 
1434
        double x,y;
 
1435
        double sxp;
 
1436
        int i;
 
1437
 
 
1438
        l->min=HUGE_VAL;
 
1439
        l->max=-HUGE_VAL;
 
1440
        l->s1=0;
 
1441
        l->sx=0;
 
1442
        l->sy=0;
 
1443
        l->sxy=0;
 
1444
        l->sxx=0;
 
1445
        for(i=0;i<l->n;i++){
 
1446
                x=l->x0+i*l->dx;
 
1447
                y=l->y_data[i];
 
1448
 
 
1449
                if(isnan(y))continue;
 
1450
 
 
1451
                if(l->y_data[i]<l->min)l->min=l->y_data[i];
 
1452
                if(l->y_data[i]>l->max)l->max=l->y_data[i];
 
1453
                l->s1+=1;
 
1454
                l->sx+=x;
 
1455
                l->sy+=y;
 
1456
                l->sxy+=x*y;
 
1457
                l->sxx+=x*x;
 
1458
        }
 
1459
        sxp=l->sxx-l->sx*l->sx/l->s1;
 
1460
        
 
1461
        l->ave_x=l->sx/l->s1;
 
1462
        l->ave_y=l->sy/l->s1;
 
1463
        l->slope=(l->s1*l->sxy-l->sx*l->sy)/(l->s1*l->sxx-l->sx*l->sx);
 
1464
        l->err_slope=l->yerr/sqrt(sxp);
 
1465
        l->err_ave_y=l->yerr/sqrt(l->s1);
 
1466
 
 
1467
        calculate_residuals(l);
 
1468
 
 
1469
        return 0;
 
1470
}
 
1471
 
 
1472
int calculate_residuals(linear_fit_t *l)
 
1473
{
 
1474
        double x,y;
 
1475
        double res,sum_res2;
 
1476
        int i;
 
1477
 
 
1478
        sum_res2=0;
 
1479
        for(i=0;i<l->n;i++){
 
1480
                x=l->x0+i*l->dx-l->ave_x;
 
1481
                y=l->y_data[i];
 
1482
 
 
1483
                if(isnan(y))continue;
 
1484
 
 
1485
                res=l->ave_y+l->slope*x-y;
 
1486
 
 
1487
                sum_res2+=res*res;
 
1488
        }
 
1489
        l->S_min=sum_res2/(l->yerr*l->yerr);
 
1490
        l->dof=l->s1-2;
 
1491
 
 
1492
        return 0;
 
1493
}
 
1494
 
 
1495
double linear_fit_func_y(linear_fit_t *l,double x)
 
1496
{
 
1497
        return l->ave_y+l->slope*(x-l->ave_x);
 
1498
}
 
1499
 
 
1500
double linear_fit_func_x(linear_fit_t *l,double y)
 
1501
{
 
1502
        return l->ave_x+(y-l->ave_y)/l->slope;
 
1503
}
 
1504
 
 
1505
void dump_curve(linear_fit_t *l)
 
1506
{
 
1507
        static int dump_number=0;
 
1508
        double x,y;
 
1509
        int i;
 
1510
 
 
1511
        printf("start dump %d\n",dump_number);
 
1512
        for(i=0;i<l->n;i++){
 
1513
                x=l->x0+i*l->dx-l->ave_x;
 
1514
                y=l->y_data[i];
 
1515
                printf("D%d: %d %g %g %g\n",dump_number,i,y,
 
1516
                        l->ave_y+l->slope*x,
 
1517
                        l->ave_y+l->slope*x-y);
 
1518
        }
 
1519
        printf("end dump\n");
 
1520
        dump_number++;
 
1521
}
 
1522
 
 
1523
 
 
1524
/* printing of scientific numbers (with errors) */
 
1525
 
 
1526
int sci_sprint(char *s,double x,double y)
 
1527
{
 
1528
        int errsig;
 
1529
        int maxsig;
 
1530
        int sigfigs;
 
1531
        double mantissa;
 
1532
        double error;
 
1533
        double mindigit;
 
1534
 
 
1535
        errsig = floor(log10(y));
 
1536
        maxsig = floor(log10(x));
 
1537
        mindigit = pow(10,errsig);
 
1538
 
 
1539
        if(maxsig<errsig)maxsig=errsig;
 
1540
 
 
1541
        sigfigs = maxsig-errsig+2;
 
1542
 
 
1543
        mantissa = x*pow(10,-maxsig);
 
1544
        error = y*pow(10,-errsig+1);
 
1545
 
 
1546
        return sprintf(s,"%0.*f(%2.0f)e%d",sigfigs-1,mantissa,error,maxsig);
 
1547
}
 
1548
 
 
1549
int sci_sprint_alt(char *s,double x,double y)
 
1550
{
 
1551
        int errsig;
 
1552
        int maxsig;
 
1553
        int sigfigs;
 
1554
        double mantissa;
 
1555
        double error;
 
1556
        double mindigit;
 
1557
 
 
1558
        errsig = floor(log10(fabs(y)));
 
1559
        maxsig = floor(log10(fabs(x)));
 
1560
        mindigit = pow(10,errsig);
 
1561
 
 
1562
        if(maxsig<errsig)maxsig=errsig;
 
1563
 
 
1564
        sigfigs = maxsig-errsig+2;
 
1565
 
 
1566
        mantissa = x*pow(10,-maxsig);
 
1567
        error = y*pow(10,-errsig+1);
 
1568
 
 
1569
        if(isnan(x)){
 
1570
                return sprintf(s,"%g",x);
 
1571
        }
 
1572
        if(errsig==1 && maxsig<4 && maxsig>1){
 
1573
                return sprintf(s,"%0.0f(%2.0f)",x,error);
 
1574
        }
 
1575
        if(maxsig<=0 && maxsig>=-2){
 
1576
                return sprintf(s,"%0.*f(%2.0f)",sigfigs-1-maxsig,
 
1577
                        mantissa*pow(10,maxsig),error);
 
1578
        }
 
1579
        return sprintf(s,"%0.*f(%2.0f)e%d",sigfigs-1,mantissa,error,maxsig);
 
1580
}
 
1581
 
 
1582
double very_low_target( comedi_t *dev, unsigned int subdevice,
 
1583
        unsigned int channel, unsigned int range )
 
1584
{
 
1585
        comedi_range *range_ptr;
 
1586
        int max_data;
 
1587
 
 
1588
        range_ptr = comedi_get_range( dev, subdevice, channel, range );
 
1589
        assert( range_ptr != NULL );
 
1590
        max_data = comedi_get_maxdata( dev, subdevice, 0 );
 
1591
        assert( max_data > 0 );
 
1592
 
 
1593
        return comedi_to_phys( 1, range_ptr, max_data ) / 2.0;
 
1594
}
 
1595
 
 
1596
double fractional_offset( calibration_setup_t *setup, int subdevice,
 
1597
        unsigned int channel, unsigned int range, int obs )
 
1598
{
 
1599
        comedi_range *range_ptr;
 
1600
        double target;
 
1601
        double reading;
 
1602
        unsigned int chanspec;
 
1603
        new_sv_t sv;
 
1604
 
 
1605
        if( subdevice < 0 || obs < 0 ) return 0.0;
 
1606
 
 
1607
        chanspec = setup->observables[obs].observe_insn.chanspec;
 
1608
        target = setup->observables[obs].target;
 
1609
 
 
1610
        range_ptr = comedi_get_range( setup->dev, subdevice, channel, range );
 
1611
        assert( range_ptr != NULL );
 
1612
 
 
1613
        comedi_set_global_oor_behavior( COMEDI_OOR_NUMBER );
 
1614
        preobserve( setup, obs);
 
1615
 
 
1616
        my_sv_init( &sv, setup, setup->ad_subdev, chanspec );
 
1617
        new_sv_measure( setup->dev, &sv );
 
1618
        reading = sv.average;
 
1619
 
 
1620
        return ( reading - target ) / ( range_ptr->max - range_ptr->min );
 
1621
}
 
1622
 
 
1623
double get_tolerance( calibration_setup_t *setup, int subdevice,
 
1624
        double num_bits )
 
1625
{
 
1626
        int maxdata;
 
1627
 
 
1628
        if( subdevice < 0 ) return INFINITY;
 
1629
 
 
1630
        maxdata = comedi_get_maxdata( setup->dev, subdevice, 0 );
 
1631
        assert( maxdata > 0 );
 
1632
        return num_bits / maxdata;
 
1633
}