summaryrefslogtreecommitdiff
path: root/utils/raspberrypi/ctt/ctt_geq.py
blob: 2aa668f133a80371929a83a9aa91268d7955c279 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
# SPDX-License-Identifier: BSD-2-Clause
#
# Copyright (C) 2019, Raspberry Pi (Trading) Limited
#
# ctt_geq.py - camera tuning tool for GEQ (green equalisation)

from ctt_tools import *
import matplotlib.pyplot as plt
import scipy.optimize as optimize


"""
Uses green differences in macbeth patches to fit green equalisation threshold
model. Ideally, all macbeth chart centres would fall below the threshold as
these should be corrected by geq.
"""
def geq_fit(Cam, plot):
    imgs = Cam.imgs
    """
    green equalisation to mitigate mazing.
    Fits geq model by looking at difference
    between greens in macbeth patches
    """
    geqs = np.array([geq(Cam, Img)*Img.againQ8_norm for Img in imgs])
    Cam.log += '\nProcessed all images'
    geqs = geqs.reshape((-1, 2))
    """
    data is sorted by green difference and top half is selected since higher
    green difference data define the decision boundary.
    """
    geqs = np.array(sorted(geqs, key=lambda r: np.abs((r[1]-r[0])/r[0])))

    length = len(geqs)
    g0 = geqs[length//2:, 0]
    g1 = geqs[length//2:, 1]
    gdiff = np.abs(g0-g1)
    """
    find linear fit by minimising asymmetric least square errors
    in order to cover most of the macbeth images.
    the philosophy here is that every macbeth patch should fall within the
    threshold, hence the upper bound approach
    """
    def f(params):
        m, c = params
        a = gdiff - (m*g0+c)
        """
        asymmetric square error returns:
            1.95 * a**2 if a is positive
            0.05 * a**2 if a is negative
        """
        return(np.sum(a**2+0.95*np.abs(a)*a))

    initial_guess = [0.01, 500]
    """
    Nelder-Mead is usually not the most desirable optimisation method
    but has been chosen here due to its robustness to undifferentiability
    (is that a word?)
    """
    result = optimize.minimize(f, initial_guess, method='Nelder-Mead')
    """
    need to check if the fit worked correectly
    """
    if result.success:
        slope, offset = result.x
        Cam.log += '\nFit result: slope = {:.5f} '.format(slope)
        Cam.log += 'offset = {}'.format(int(offset))
        """
        optional plotting code
        """
        if plot:
            x = np.linspace(max(g0)*1.1, 100)
            y = slope*x + offset
            plt.title('GEQ Asymmetric \'Upper Bound\' Fit')
            plt.plot(x, y, color='red', ls='--', label='fit')
            plt.scatter(g0, gdiff, color='b', label='data')
            plt.ylabel('Difference in green channels')
            plt.xlabel('Green value')

        """
        This upper bound asymmetric gives correct order of magnitude values.
        The pipeline approximates a 1st derivative of a gaussian with some
        linear piecewise functions, introducing arbitrary cutoffs. For
        pessimistic geq, the model parameters have been increased by a
        scaling factor/constant.

        Feel free to tune these or edit the json files directly if you
        belive there are still mazing effects left (threshold too low) or if you
        think it is being overcorrected (threshold too high).
        We have gone for a one size fits most approach that will produce
        acceptable results in most applications.
        """
        slope *= 1.5
        offset += 201
        Cam.log += '\nFit after correction factors: slope = {:.5f}'.format(slope)
        Cam.log += ' offset = {}'.format(int(offset))
        """
        clamp offset at 0 due to pipeline considerations
        """
        if offset < 0:
            Cam.log += '\nOffset raised to 0'
            offset = 0
        """
        optional plotting code
        """
        if plot:
            y2 = slope*x + offset
            plt.plot(x, y2, color='green', ls='--', label='scaled fit')
            plt.grid()
            plt.legend()
            plt.show()

        """
    the case where for some reason the fit didn't work correctly

    Transpose data and then least squares linear fit. Transposing data
    makes it robust to many patches where green difference is the same
    since they only contribute to one error minimisation, instead of dragging
    the entire linear fit down.
    """

    else:
        print('\nError! Couldn\'t fit asymmetric lest squares')
        print(result.message)
        Cam.log += '\nWARNING: Asymmetric least squares fit failed! '
        Cam.log += 'Standard fit used could possibly lead to worse results'
        fit = np.polyfit(gdiff, g0, 1)
        offset, slope = -fit[1]/fit[0], 1/fit[0]
        Cam.log += '\nFit result: slope = {:.5f} '.format(slope)
        Cam.log += 'offset = {}'.format(int(offset))
        """
        optional plotting code
        """
        if plot:
            x = np.linspace(max(g0)*1.1, 100)
            y = slope*x + offset
            plt.title('GEQ Linear Fit')
            plt.plot(x, y, color='red', ls='--', label='fit')
            plt.scatter(g0, gdiff, color='b', label='data')
            plt.ylabel('Difference in green channels')
            plt.xlabel('Green value')
        """
        Scaling factors (see previous justification)
        The model here will not be an upper bound so scaling factors have
        been increased.
        This method of deriving geq model parameters is extremely arbitrary
        and undesirable.
        """
        slope *= 2.5
        offset += 301
        Cam.log += '\nFit after correction factors: slope = {:.5f}'.format(slope)
        Cam.log += ' offset = {}'.format(int(offset))

        if offset < 0:
            Cam.log += '\nOffset raised to 0'
            offset = 0

        """
        optional plotting code
        """
        if plot:
            y2 = slope*x + offset
            plt.plot(x, y2, color='green', ls='--', label='scaled fit')
            plt.legend()
            plt.grid()
            plt.show()

    return round(slope, 5), int(offset)


""""
Return green channels of macbeth patches
returns g0, g1 where
> g0 is green next to red
> g1 is green next to blue
"""
def geq(Cam, Img):
    Cam.log += '\nProcessing image {}'.format(Img.name)
    patches = [Img.patches[i] for i in Img.order][1:3]
    g_patches = np.array([(np.mean(patches[0][i]), np.mean(patches[1][i])) for i in range(24)])
    Cam.log += '\n'
    return(g_patches)
strerror(ret); count = errorIdx - 1; ret = errorIdx; } updateControls(ctrls, controlInfo, v4l2Ctrls, count); return ret; } /** * \brief Write controls to the device * \param[in] ctrls The list of controls to write * * This method writes the value of all controls contained in \a ctrls, and * stores the values actually applied to the device in the corresponding * \a ctrls entry. * * If any control in \a ctrls is not supported by the device, is disabled (i.e. * has the V4L2_CTRL_FLAG_DISABLED flag set), is read-only, is a * compound control, or if any other error occurs during validation of * the requested controls, no control is written and this method returns * -EINVAL. * * If an error occurs while writing the controls, the index of the first * control that couldn't be written is returned. All controls below that index * are written and their values are updated in \a ctrls, while all other * controls are not written and their values are not changed. * * \return 0 on success or an error code otherwise * \retval -EINVAL One of the control is not supported or not accessible * \retval i The index of the control that failed */ int V4L2Device::setControls(V4L2ControlList *ctrls) { unsigned int count = ctrls->size(); if (count == 0) return 0; const V4L2ControlInfo *controlInfo[count]; struct v4l2_ext_control v4l2Ctrls[count]; memset(v4l2Ctrls, 0, sizeof(v4l2Ctrls)); for (unsigned int i = 0; i < count; ++i) { const V4L2Control *ctrl = ctrls->getByIndex(i); const auto iter = controls_.find(ctrl->id()); if (iter == controls_.end()) { LOG(V4L2, Error) << "Control '" << ctrl->id() << "' not found"; return -EINVAL; } const V4L2ControlInfo *info = &iter->second; controlInfo[i] = info; v4l2Ctrls[i].id = info->id(); /* Set the v4l2_ext_control value for the write operation. */ switch (info->type()) { case V4L2_CTRL_TYPE_INTEGER64: v4l2Ctrls[i].value64 = ctrl->value().get<int64_t>(); break; default: /* * \todo To be changed when support for string and * compound controls will be added. */ v4l2Ctrls[i].value = ctrl->value().get<int32_t>(); break; } } struct v4l2_ext_controls v4l2ExtCtrls = {}; v4l2ExtCtrls.which = V4L2_CTRL_WHICH_CUR_VAL; v4l2ExtCtrls.controls = v4l2Ctrls; v4l2ExtCtrls.count = count; int ret = ioctl(VIDIOC_S_EXT_CTRLS, &v4l2ExtCtrls); if (ret) { unsigned int errorIdx = v4l2ExtCtrls.error_idx; /* Generic validation error. */ if (errorIdx == 0 || errorIdx >= count) { LOG(V4L2, Error) << "Unable to set controls: " << strerror(ret); return -EINVAL; } /* A specific control failed. */ LOG(V4L2, Error) << "Unable to set control " << errorIdx << ": " << strerror(ret); count = errorIdx - 1; ret = errorIdx; } updateControls(ctrls, controlInfo, v4l2Ctrls, count); return ret; } /** * \brief Perform an IOCTL system call on the device node * \param[in] request The IOCTL request code * \param[in] argp A pointer to the IOCTL argument * \return 0 on success or a negative error code otherwise */ int V4L2Device::ioctl(unsigned long request, void *argp) { /* * Printing out an error message is usually better performed * in the caller, which can provide more context. */ if (::ioctl(fd_, request, argp) < 0) return -errno; return 0; } /** * \fn V4L2Device::deviceNode() * \brief Retrieve the device node path * \return The device node path */ /** * \fn V4L2Device::fd() * \brief Retrieve the V4L2 device file descriptor * \return The V4L2 device file descriptor, -1 if the device node is not open */ /* * \brief List and store information about all controls supported by the * V4L2 device */ void V4L2Device::listControls() { struct v4l2_query_ext_ctrl ctrl = {}; /* \todo Add support for menu and compound controls. */ while (1) { ctrl.id |= V4L2_CTRL_FLAG_NEXT_CTRL; if (ioctl(VIDIOC_QUERY_EXT_CTRL, &ctrl)) break; if (ctrl.type == V4L2_CTRL_TYPE_CTRL_CLASS || ctrl.flags & V4L2_CTRL_FLAG_DISABLED) continue; switch (ctrl.type) { case V4L2_CTRL_TYPE_INTEGER: case V4L2_CTRL_TYPE_BOOLEAN: case V4L2_CTRL_TYPE_MENU: case V4L2_CTRL_TYPE_BUTTON: case V4L2_CTRL_TYPE_INTEGER64: case V4L2_CTRL_TYPE_BITMASK: case V4L2_CTRL_TYPE_INTEGER_MENU: break; /* \todo Support compound controls. */ default: LOG(V4L2, Debug) << "Control type '" << ctrl.type << "' not supported"; continue; } controls_.emplace(std::piecewise_construct, std::forward_as_tuple(ctrl.id), std::forward_as_tuple(ctrl)); } } /* * \brief Update the value of the first \a count V4L2 controls in \a ctrls using * values in \a v4l2Ctrls * \param[inout] ctrls List of V4L2 controls to update * \param[in] controlInfo List of V4L2 control information * \param[in] v4l2Ctrls List of V4L2 extended controls as returned by the driver * \param[in] count The number of controls to update */ void V4L2Device::updateControls(V4L2ControlList *ctrls, const V4L2ControlInfo **controlInfo, const struct v4l2_ext_control *v4l2Ctrls, unsigned int count) { for (unsigned int i = 0; i < count; ++i) { const struct v4l2_ext_control *v4l2Ctrl = &v4l2Ctrls[i]; const V4L2ControlInfo *info = controlInfo[i]; V4L2Control *ctrl = ctrls->getByIndex(i); switch (info->type()) { case V4L2_CTRL_TYPE_INTEGER64: ctrl->value().set<int64_t>(v4l2Ctrl->value64); break; default: /* * \todo To be changed when support for string and * compound controls will be added. */ ctrl->value().set<int32_t>(v4l2Ctrl->value); break; } } } } /* namespace libcamera */