© ESO, 2017

1. Introduction

Solar rotation can be generally determined with three methods: the tracer method, the spectroscopic method, and the helioseismology method (Stix 2004). Sunspots and sunspot groups represent one of the most commonly used tracers in the literature. The main reason is the availability of long- term data sets from various observatories: the Greenwich photoheliographic results (GPR) data set (Balthasar & Wöhl 1980; Arévalo et al. 1982; Balthasar et al. 1986b); the Mt. Wilson data set (Howard 1984; Gilman & Howard 1984; Hathaway & Wilson 1990); the extended Greenwich results (EGR) data set (Pulkkinen & Tuominen 1998; Javaraiah 2003; Zuccarello & Zappalá 2003; Javaraiah et al. 2005; Javaraiah & Ulrich 2006; Brajša et al. 2006, 2007; Sudar et al. 2014); the Kodaikanal data set (Gupta et al. 1999); the Debrecen photoheliographic data (DPD; Sudar et al. 2017); and the Kanzelhöhe Observatory for Solar and Environmental Research (KSO) data set (Lustig 1982, 1983; Balthasar & Fangmeier 1988; Poljančić Beljan et al. 2014; Pötzi et al. 2016).

The KSO provides daily multispectral synoptic observations of the Sun using several telescopes (Veronig & Pötzi 2016). However, especially the sunspot drawings and white light images seem to be neglected by the scientific community. Previous studies performed on solar drawings from the KSO, in which sunpots and sunspot groups were analysed, involve data obtained through the year 1985 (Lustig 1982, 1983; Lustig & Dvorak 1984; Balthasar et al. 1986a; Hanslmeier & Lustig 1986; Lustig & Hanslmeier 1987; Balthasar & Fangmeier 1988; Lustig & Wöhl 1991). There are a few papers that use more recent sunspot drawings (Temmer et al. 2006; Poljančić et al. 2010; Poljančić et al. 2011; Poljančić Beljan et al. 2014). Temmer et al. (2006) provide the catalog of the hemispheric sunspot numbers and investigate the north – south activity asymmetries, but do not deal with the solar rotation.

Poljančić et al. (2010) checked the precision ofthe Solar Optical Observing Network/United States Air Force/National Oceanic and Atmospheric Administration (SOON/USAF/ NOAA) data set by comparing it with the GPR data set. The KSO data set (heliographic coordinates determined from the sunspot drawings for the years 1972 and 1993) were used as a reference for the comparison of sunspot position measurements. The SOON/USAF/NOAA data were found to be somewhat less accurate than the GPR data. Poljančić et al. (2011) expanded that analysis on several more data sets, again by using the KSO data set as a reference, for the comparison of sunspot position measurements. In almost all cases the latitude, longitude, and synodic angular velocity differences between the KSO and other data sets are the smallest for DPD, GPR, and SOON/USAF/NOAA. Since the GPR and SOON/USAF/NOAA data sets have already been widely used for the analysis of differential rotation of sunspot groups and its temporal variation, we decided that the next step will be to process the DPD and KSO data sets. An analysis concerning the KSO data set is presented in this work, while an analysis of the DPD data set is made in Sudar et al. (2017). Poljančić Beljan et al. (2014) calculated differential rotation parameters and showed differential rotation profiles for solar cycles Nos. 20 and 22 for the KSO data in a preliminary form.

Our analysis is a continuation of the investigation done by Lustig (1983). We processed the KSO sunspot drawings for solar cycles Nos. 20–23 (1964–2008), and the KSO white light images for solar cycle No. 24 (2009–2016). For the first time, the whole solar cycle No. 21 is examined using the KSO data (previously, it was just partially processed), as well as solar cycles Nos. 22–24. Solar cycle No. 20 is where our work and the work of Lustig (1983) overlap and is therefore suitable for comparison of the results.

In this paper we study the solar differential rotation and we make a comparison of the differential rotation parameters collected from different sources. In addition, we discuss the north-south asymmetry of the solar rotation profile.

2. Instrumentation and measurements

Fig. 1 Projection device for making the draft of sunspot drawings. The mirrors in the light path of the refractor have reversed the drawing, i.e., west is left. The enlargement of the solar disk to a diameter of about 25 cm corresponds to a magnification factor of about 100.

The KSO belongs to the University of Graz, Austria¹. Sunspot drawings are made with a refractor (d/f = 110/1650 mm) using the projection system that enlarges the disk image to 25 cm in diameter (see Fig. 1). The drawings are updated online almost every day, depending on weather conditions. In general, weather conditions are good and continuous operation yields about 300 days of observations per year. The KSO archive consists of more than 60 yr of sunspot drawings (Pötzi et al. 2016). Otruba (2006) described the solar monitoring program at Kanzelhöhe Observatory and the sunspot drawings in detail.

The white light images are made by a refractor (d/f = 130/1950 mm) using a filter for 546 nm with a band pass of 10 nm. During the period 1989–2007 white light images were recorded on photographic film. Usually, three white light images per day were produced. Until recently, these images were accessible only on transparent sheet films, but the scanning process is proceeding and they are now available online for the years after 1992 (see Fig. 2). The films are scanned with a photo scanner for transparent film material and as an output FITS and JPEG files are produced (Pötzi 2010).

In July 2007 the system was replaced by the KSO Photosphere Digital Camera (KPDC), Pulnix TM-4100CL (Pötzi 2010). To improve the image quality a 12 bit 2048 × 2048 pixel Pulnix RM-4200GE camera was installed in August 2015.

Fig. 2 Scanned white light image with added NOAA numbers and a grid that shows P and B0 angles. The image is not rotated, but an interface for de-rotating the images into the east-west direction is available (2002-08-19 06:10:41 UT).

3. Determination of the heliographic positions

Our investigation covers the time period 1964–April 2016. We used two procedures to determine the heliographic positions of sunspot groups: an interactive procedure, using KSO sunspot drawings for the time period 1964–2008, and an automatic procedure, using KSO white light images for the time period 2009–2016.

In order to avoid solar limb effects leading to high position uncertainties, we limited the data to ±58 deg in central meridian distance (CMD) which covers about 85% of the projected solar radius (Balthasar et al. 1986b). With this cutoff we obtained a sample of 12 152 sunspot groups that correspond to some 70 000 individual positions of sunspot groups. Recurrent sunspot groups are counted as many times as they appear.

3.1. Interactive procedure (1964–2008, solar cycles Nos. 20–23)

A software package called Sungrabber² (Hržina et al. 2007) that determines the position of tracers is used for the interactive procedure. The area weighted centers of sunspot groups are estimated by the naked eye, thus giving more importance (by moving the center of gravity) to the umbra and penumbra that are more pronounced, i.e., occupy a larger area.

In order to find out if different observers could affect the measurements, we performed a test: two observers independently determined the positions of 10 single H- and J-type sunspot groups (see Poljančić et al. 2010, Fig. 3). Differences between measurements of the two observers are negligible (on average ≈0.1 deg, always less then 0.2 deg) when compared with the differences between various observatories (≈0.5 deg). Based on these results all interactive coordinate determinations in the present work were done by two different observers using the Sungrabber software.

Sunspot groups were identified with the help of the GPR and DPD databases (Baranyi et al. 2016; Győri et al. 2017). The software Sungrabber offers the possibility of loading the GPR or DPD sunspot group heliographic positions, as well as Greenwich or NOAA/USAF sunspot group numbers, and mark them on the sunspot drawing used.

3.2. Automatic procedure (2009–2016, solar cycle No. 24)

Based on the algorithm by Watson & Fletcher (2011) sunspot groups and their properties (size, umbra, penumbra, position) are identified by morphological image processing of KSO white light images. This data is prepared every observing day by KSO and available over KSO ftp server³ as raw fits images⁴ and data files⁵. Based on the assumption that images of better quality show more detail, we selected for each single day the fits file with the corresponding data file that contained the largest amount of information.

The data files provide the sunspot group center of gravity in pixel coordinates. Only the umbra pixels are used to calculate the center of gravity. For each sunspot group, a few separate umbra/penumbra pixel coordinates are also available, determined as the mean position of the umbra/penumbra pixels, which are only area weighted.

From the fits header P and B₀ were extracted, and from the corresponding data file solar radii and X and Y pixel coordinates were extracted. Using the calculations from Meeus (1991), we calculated heliographic coordinates of sunspot groups.

We also tried to apply the automatic method on white light images recorded on photographic film material (1989–2007), but the results were not very reliable. The calculated equatorial rotation velocities ranged from less than 14 deg/day to more than 15 deg/day. The reason for these uncertainties is the strong variation of the density (blackening) due to the developing procedure of the film material. These variations lead to some very faint images where the solar limb is very difficult to detect, leading to incorrect solar radii and as a consequence to incorrect sunspot positions.

4. Determining the rotation velocities

In order to calculate the synodic angular rotation velocities we used two methods: the daily-shift method (DS), where the synodic rotation velocities were calculated from the daily differences of the CMD and the elapsed time t${\mathit{\omega }}_{\mathrm{syn}}\mathrm{=}\frac{\mathrm{\Delta }\mathit{CMD}}{\mathrm{\Delta }\mathit{t}}\mathit{,}$(1)and the robust linear least-squares fit method (rLSQ), where the synodic rotation velocities were calculated by fitting a line to the measured positions in time CMD(t) for each tracer. The synodic rotation velocity corresponds to the slope of the fit. Since the measured data occasionally have outliers, due to false identification and other reasons, we used a robust fit, namely iteratively reweighted least-squares with Huber’s t weighting function (Huber 1981). When applying the rLSQ method, the rotation velocities were calculated by fitting a line to at least three data points. If the number of position measurements, i.e., data points, is calculated for each sunspot group, the median value is 5, while the maximum value reaches 11. We note that even in the cases where only three data points is used by the rLSQ method, it is still more than the two data points used when applying the DS method.

Angular rotation velocities were converted from synodic to sidereal using the procedure described by Roša et al. (1995) and by Brajša et al. (2002a), which was improved by Skokić et al. (2014). The calculated sidereal velocities, ω, were used in the least-squares fitting to the solar differential rotation law $\mathit{\omega }\mathrm{\left(}\mathit{b}\mathrm{\right)}\mathrm{=}\mathit{A}\mathrm{+}\mathit{B}{\sin }^{\mathrm{2}}\mathit{b,}$(2)where b is the heliographic latitude, and A and B are the solar differential rotation parameters. In both methods, DS and rLSQ, we assigned the velocity to the latitude and time of the first measurement of position. When average latitudes are used, the false meridional flows can be detected (Olemskoy & Kitchatinov 2005), but for rotation analysis it probably represents a negligible effect.

5. Results

5.1. Comparison of automatic and interactive procedures

A direct comparison of the interactive and the automatic procedure can only be presented for the year 2014, the only year for which both procedures were applied. This year was selected as it represents the maximum of solar cycle No. 24 after a deep minimum, when almost no sunspots and thus data points were available, and the first maximum covered by the KSO Photosphere Digital Camera. We note that the years near the solar minimum usually lead to large statistical errors as the data set becomes very small, e.g., the year 2009 with 260 spotless days.

The differential rotation parameters A and B determined for the year 2014 with both procedures and the velocities obtained with DS and rLSQ methods for both hemispheres together are presented in Table 1. Corresponding differential rotation profiles are shown in Fig. 3. More details are given in Sect. 6.1.

Fig. 3 Differential rotation profiles for the year 2014 (see Table 1). Sidereal rotation velocity is denoted by ω and heliographic latitude by b.

5.2. Differential rotation from the KSO data

Using both procedures, interactive and automatic, five solar cycles (20–24) from 1964 until 2016 were used to obtain the differential rotation parameters. Both solar hemispheres were treated together and separately.

With only a CMD cutoff we obtained a sample of 41 125 (33 817 via DS and 7308 via rLSQ) calculated sidereal rotation velocities. Rotational velocity outliers resulting from misidentification of sunspot groups in subsequent images can be filtered out by applying the velocity filter 8–19 deg/day (Brajša et al. 2002b; Vršnak et al. 2003; Sudar et al. 2014, 2015). After selecting only sidereal rotation velocities higher than 8 deg/day and lower than 19 deg/day, the number of calculated sidereal rotation velocities was reduced to 40 480 (33 213 via DS and 7267 via rLSQ). This means that 1.82% of the calculated velocities were eliminated using the DS method and 0.56% using the rLSQ method. This also means that 51% of the reduced calculated sidereal rotation velocities belong to the northern hemisphere, while 49% belong to the southern hemisphere (valid for both the DS and the rLSQ).

In Table 2 we present the results of fittings for solar cycles Nos. 20–24 separately (DS: rows 1–5 for N+S, rows 8–12 for N, rows 15–19 for S; rLSQ: rows 22–26 for N+S, rows 29–33 for N, rows 36–40 for S); the results for the part of the data set using the interactive procedure covering solar cycles Nos. 20–23 (DS: rows 6, 13, 20; rLSQ: rows 27, 34, 41) and the results for the whole data set using the interactive and automatic method together covering solar cycles Nos. 20–24 (DS: rows 7, 14, 21; rLSQ: rows 28, 35, 42). All results are calculated using both methods (DS and rLSQ), for both hemispheres together (N+S), and for the northern hemisphere (N) and the southern hemisphere (S) separately. The years specified in Table 2 (time period) represent the starting and ending epochs of the corresponding solar cycles. They are taken from Brajša et al. (2009), except the beginning of solar cycle No. 24, which was taken from SILSO World Data Center (2015), Royal Observatory of Belgium, Brussels⁶.

5.3. Comparison with other data

In this section we present a comparison of the differential rotation parameters A and B collected from different sources, mostly those already mentioned in Sect. 1 as interesting (GPR, EGR, DPD, KSO) dealing only with sunspots and sunspot groups as tracers (Table 3). Table 3 consists of three parts: the upper part (rows 1–21) gives the results for both solar hemispheres together, the middle part (rows 22–32) gives the results for the northern hemisphere, and the lower part (rows 33–43) gives the results for the southern hemisphere. Figure 4 only shows the differential rotation profiles for the results derived from the KSO data set (Table 3, rows 17–21). Figure 5 shows the differential rotation profiles for several data sets and both hemispheres together (see corresponding rows in Table 3), using only sunspot groups (Note (a)) or sunspots and sunspot groups (Note (c)) as tracers. The differential rotation profiles derived from the KSO data set Lustig (1983; Fig. 5, red line) and the Mt. Wilson data set Howard et al. (1984; Fig. 5, pink dashed line) show significantly lower values when compared to other results. More details are given in Sect. 6.2.

Fig. 4 Differential rotation profiles calculated by different authors for the KSO data set and both hemispheres together (see corresponding rows in Table 3). Sidereal rotation velocity is denoted by ω and heliographic latitude by b.

5.4. North–south asymmetry of the solar rotation

The differential rotation parameters derived in the present paper for the whole time period (1964–2016) and for the northern and southern hemispheres separately are given in Table 2 (DS method – rows 14 and 21; rLSQ method – rows 35 and 42). Corresponding differential rotation profiles are shown in Fig. 6. More details are given in Sect. 6.3.

6. Discussion

6.1. Automatic vs. interactive procedure, DS vs. rLSQ method

If we compare the results for the year 2014 determined by the interactive and automatic procedures and by the DS method (Table 1, rows 1 and 2), for both differential rotation parameters the results are within the 1 common σ, which also holds for the rLSQ results (Table 1, rows 3 and 4). However, both procedures, the interactive and the automatic, show large standard errors, especially for parameter B. A comparison of the standard errors in Table 1 (rows 1 and 2, rows 3 and 4) between the two procedures reveals that the automatic procedure is less accurate but can process more data. The application of the interactive procedure is rather time consuming and the results presented here are a combination of the efforts of two co-authors over almost two years of measurements.

It is obvious that the differential rotation profiles for the rLSQ method coincide better (Fig. 3, black and blue lines) because the equatorial rotation velocities are almost the same and the gradients of the solar rotation do not differ as in the DS case. In the case of the automatic procedure, the rLSQ method for calculating rotational velocities (Fig. 3, blue line) is much more reliable than the DS method (Fig. 3, red line), which can be clearly seen by a comparison of the corresponding differential rotation parameters in Table 1 (rows 2 and 4). For example, parameter B calculated via the DS method yields B = (−0.98 ± 1.11) deg/day, while for the rLSQ method B = (−3.40 ± 1.32) deg/day. If relative standard errors are compared, the rLSQ method yields a value that is three times smaller than the corresponding value for the DS method.

However, if we look at the results obtained by applying only the interactive procedure (Table 2, rows 6, 13, 20, 27, 34, 41) and compare them with the results of the interactive and automatic procedures together covering the whole time period (Table 2, rows 7, 14, 21, 28, 35, 42), we see that the inclusion of the results calculated by the automatic procedure, for solar cycle No. 24, does not affect the “whole time period” result. In other words, the int (all) and int + aut (all) results from the corresponding rows in Table 2, which are listed above, are almost identical in all cases for both the DS and rLSQ methods and for both hemispheres. Thus, when the rotation velocities calculated by the automatic procedure are combined with the rotation velocities calculated by the interactive procedure for longer periods of time, the effect of the lower accuracy of the automatic procedure has no statistically significant influence on the final best fit differential rotation profile. Hence, we combine both procedures if we analyse the whole period from 1964 to 2016. No separation between an interactive procedure before 2009 and an automatic procedure after 2009 has to be made.

For white light images recorded on photographic film material before 2007, which are scanned, the automatic procedure did not yield reliable results, and we were forced to process this data using the interactive procedure and sunspot drawings. So, for the years before 2007, it is possible to apply the automatic procedure after elimination of the possible errors that appear during the scanning procedure. However, for the years after 2007, the automatic procedure represents a useful tool for processing the white light images.

If we look at Table 2 and compare the DS and rLSQ measurements derived for the whole periods (both hemispheres – rows 7 and 28; northern hemisphere – rows 14 and 35; southern hemisphere – rows 21 and 42), the errors for both differential rotation parameters in all cases (northern, southern, or both hemispheres) are within the 1 common σ. In some cases there are very small (i.e., almost negligible) differences. This means that the results of the two different methods used to determine the rotation velocities for the whole time period (daily shift – DS and robust linear least-squares fit – rLSQ) can be considered almost identical.

Fig. 5 Differential rotation profiles calculated by different authors for several data sets (GPR, EGR, Mt. Wilson, DPD, KSO) and both hemispheres together (see corresponding rows in Table 3). Sidereal rotation velocity is denoted by ω and heliographic latitude by b.

6.2. Comparison with other data sets

In Poljančić et al. (2011) a comparison of sunspot position measurements for several data sets was made, separately for single H and J and complex sunspot groups types. The mean absolute differences between the KSO and DPD positions are lower in all cases – except for the longitude difference for complex sunspot groups – than the mean absolute difference between the KSO and GPR positions. As the KSO and GPR comparison was made for the year 1972, it may be possible that the KSO data is influenced by the erroneous solar image radii (Balthasar et al. 1984), thus causing a larger difference than that between the KSO and the DPD. The KSO and DPD comparison was made for the year 1993, and was no longer influenced by the previously mentioned unwanted effect, because the KSO observing and reduction procedures were improved after the 1980s (Balthasar et al. 1984).

We compared the differential rotation parameters derived in the present work from KSO for both the DS and the rLSQ methods (Table 3, rows 20 and 21) with values for the EGR (Sudar et al. 2014; Table 3, row 8) and values for the DPD (Sudar et al. 2017; Table 3, row 16), which were also obtained by tracing sunspot groups and covering almost the same time period. Values for the differential rotation parameter A are within the 1 common σ, while values for the differential rotation parameter B match within the 2 common σ. It is well known that the GPR represents a homogeneous data set with high accuracy (Balthasar & Wöhl 1980; Arévalo et al. 1982; Willis et al. 2013a,b; Erwin et al. 2013), as well as the DPD as its continuation (Sudar et al. 2017). Therefore, we feel free to conclude that the KSO data used in this paper (from 1964 till present) is of comparable accuracy to the DPD or GPR data sets, and suitable for the investigation of the solar rotation, although there is some indication of lower accuracy in the past, especially before the 1960s (Balthasar et al. 1984).

Except for the already mentioned overlapping with Sudar et al. (2014, 2017) rotation results (Table 3, rows 8 and 16), in most of the other studies listed in Table 3 the differential rotation parameters are in agreement with our results Table 3, rows 20 and 21); the exceptions are the GPR result for reccurent sunspot groups (Brajša et al. 2002a; Table 3, row 4); the Mt. Wilson result concerning only sunspot groups as tracers (Howard et al. 1984; Table 3, row 15); the Abastumani result (Khutsishvili et al. 2002; Table 3, row 13); and the KSO results (Lustig 1982, 1983; Balthasar & Fangmeier 1988; Table 3, rows 17–19).

According to Ruždjak et al. (2004, 2005) recurrent sunspots and recurrent sunspot groups show a slower rotation, which may be the reason for the low rotation values of the GPR given by Brajša et al. (2002a; Table 3, row 4). It is interesting that the lower value of the equatorial velocity A was also obtained for the Mt. Wilson data for sunspot groups (Howard et al. 1984; Table 3, row 15), but it is not clear what could be the cause. A higher value of the equatorial velocity, as well as the lower value of the differential rotation gradient was obtained with the Abastumani data (Khutsishvili et al. 2002; Table 3, row 13) compared to our analysis (Table 3, rows 20 and 21). This result, however, shows higher standard errors in comparison to other results in Table 3.

The equatorial rotation rates A derived from the KSO data (Lustig 1982, 1983; Balthasar & Fangmeier 1988; Table 3, rows 17–19) for the years before the 1980s are significantly lower than our analysis (Table 3 – rows 20 and 21) by about 0.15 deg. A large difference between differential rotation parameters B is also noticed. Balthasar & Fangmeier (1988) analysed the GPR data for the same time period (Table 3, row 7) and also noted the disagreement with the KSO data, confirming that the KSO data before the 1980s are affected by some systematic errors, already examined in Balthasar et al. (1984).

Since our investigations were conceived as a continuation of the investigations done in Lustig (1983), it is important to compare them also for the overlapping period of time (1964–1976, solar cycle No. 20). Lustig (1983) used the Stonyhurst disk overlaying technique to determine the heliographic coordinates of sunspot groups by superimposing the disk on sunspot drawing. For solar cycle No. 20, Lustig (1983) provides the solar rotation parameters only for the northern and southern hemispheres. They are listed in Table 3 (rows 30 and 41). The corresponding solar rotation parameters derived in the present paper are listed in Table 2 (rows 8 and 15 for DS, rows 29 and 36 for rLSQ). All these rotation profiles are shown in Fig. 7. The results for both differential rotation parameters and both methods, the DS and the rLSQ, are within the 1 common σ. Only the differential rotation parameter A for the DS method and southern hemisphere falls within the 2 common σ with the corresponding result from Lustig (1983). This comparison is evidence that the reliability of the KSO data set has been improved, not only for the period after the 1980s, but indeed from 1960 onwards.

Fig. 6 North – south asymmetry for the whole time period 1964–2016 and for DS and rLSQ methods (see also the corresponding rows in Table 2). Sidereal rotation velocity is denoted by ω and heliographic latitude by b. N – northern hemisphere, S – southern hemisphere.

6.3. North–south asymmetry of the solar rotation

If we compare the DS measurements derived for the whole time period (1964–2016) and for the northern and the southern hemispheres separately (Table 2, rows 14 and 21), we see that the measurements for differential rotation parameter A are within the 1 common σ, i.e., the difference is very small and statistically insignificant. Differential rotation parameter B shows the difference between the two solar hemispheres, statistically significant on 2 common σ level (Table 2, rows 14 and 21). Concerning the rLSQ measurements (Table 2, rows 35 and 42), we see that the measurements for both differential rotation parameters are within the 1 common σ.

In Lustig (1983), statistically insignificant differences were discovered for equatorial velocities (the parameter A) between the northern and southern hemispheres, while the difference in the gradients of the differential rotation (the parameter B) is 0.36 deg/day. It is interesting that in present work for the DS method this difference remains almost the same (B_N−B_S = 0.31 deg/day).

For both methods, the DS and the rLSQ, the southern solar hemisphere rotates a little bit faster, which can also be clearly seen in Fig. 6. Therefore, a barely noticeable north–south asymmetry is observed for the whole time period 1964–2016 in the present work.

7. Conclusions

We have determined the solar differential rotation by tracing sunspot groups during the period 1964–2016. We have used two procedures to determine the heliographic positions: an interactive one on the KSO sunspot drawings (1964–2008, solar cycles Nos. 20–23) and an automatic one on the KSO white light images (2009–2016, solar cycle No. 24). For the first time, the whole solar cycle No. 21 was investigated with the KSO data, as well as solar cycles Nos. 22–24. Applying a CMD cutoff at ±58 deg we obtained a sample of 12 152 sunspot groups which correspond to approximately 70 000 individual sunspot positions. The synodic angular rotation velocities were determined using two different methods, the DS and rLSQ, and then converted to the sidereal velocities. A sample of 41 125 calculated sidereal rotation velocities was obtained. Calculated sidereal velocity values have been used in the least-squares fitting to the solar differential rotation law.

The interactive and automatic procedures comparison was performed for the year 2014. It yielded 1 common σ coincidence for both differential rotation parameters and both the DS and rLSQ methods. However, in both procedures, interactive and automatic, there are large standard errors for the differential rotation parameters, especially for the differential rotation parameter B. If the comparison were performed for the entire cycle, it would yield better results as the number of calculated velocities and latitudinal coverage is higher in that case. For the automatic procedure the standard errors are even larger, indicating lower accuracy of the automatic procedure in comparison to the interactive one. Also, in the case of automatic procedure, the rLSQ method for calculating rotational velocities is much more reliable than the DS method. For the test data from 2014, the rLSQ method yields a relative standard error for the differential rotation parameter B that is three times smaller than that of the DS method. But, when the rotation velocities calculated by the automatic procedure are combined with the rotation velocities calculated by the interactive procedure for longer periods of time, the effect of the lower accuracy of the automatic procedure has no statistically significant influence on the final best fit differential rotation profile. This allows us to discuss interactive procedure before 2009 and automatic procedure after 2009 together as if only one procedure has been applied.

Fig. 7 Comparison between the present paper and Lustig (1983) differential rotation profiles calculated for solar cycle No. 20 (1964–1976) and for the northern and southern hemisphere separately, see corresponding rows in Tables 2 and 3. The sidereal rotation velocity is denoted by ω and heliographic latitude by b.

The two different methods used for the determination of the rotation velocities (daily shift – DS and robust linear least-squares fit – rLSQ) were also compared for the whole time period, and showed 1 common σ coincidence for both differential rotation parameters.

Inspection for the north – south asymmetry showed that for both methods, the DS and the rLSQ, the southern solar hemisphere rotates a little bit faster. Thereby, a barely noticeable north – south asymmetry is observed for the whole time period 1964–2016 in the present paper, very similar to previous observations for the time period 1947–1981 (Lustig 1983).

The comparison of calculated rotation profiles derived from the KSO data in this work and from all other studies collected from different sources yielded a conclusion about the KSO accuracy after the 1960s. The KSO data used in this paper for the time period from 1964 to the present is of a comparable accuracy with the DPD or GPR data sets, and suitable for the investigation of the solar rotation. We cannot guarantee that for the years prior to 1960 the results would be of the same accuracy because there are some indications about lower accuracy in the past (Balthasar et al. 1984).

The quality of the KSO sunspot drawings has gradually increased over the last 50 yr. In general, the interactive procedure of position determination is fairly accurate, but has the drawback of being very time consuming, so the much faster automatic procedure of position determination was developed. However, in the case of the automatic procedure, the rLSQ method for calculating rotational velocities is much more reliable than the DS method.

The main conclusion of this paper is that KSO provides a valuable data set with satisfactory accuracy. Although there were some systematic errors in the past, this data set is well suited to long-term studies and therefore still continued.

We plan to continue our analysis of the KSO data. The next step is presenting the temporal variation of the differential rotation and the relationship between the solar rotation and activity, as well as an analysis of meridional motions and horizontal Reynolds stress.

Acknowledgments

This work was supported in part by the Croatian Science Foundation under the project 6212 “Solar and Stellar Variability” and in part by the University of Rijeka under project number 13.12.1.3.03. It has also received funding from the SOLARNET project (312495, 2013–2017), which is an Integrated Infrastructure Initiative (I3) supported by FP7 Capacities Programme. I.P.B. thanks the Kanzelhöhe Solar Observatory for the hospitality during her stay at Kanzelhöhe.

References

All Tables

All Figures

	Fig. 1 Projection device for making the draft of sunspot drawings. The mirrors in the light path of the refractor have reversed the drawing, i.e., west is left. The enlargement of the solar disk to a diameter of about 25 cm corresponds to a magnification factor of about 100.
In the text

	Fig. 2 Scanned white light image with added NOAA numbers and a grid that shows P and B0 angles. The image is not rotated, but an interface for de-rotating the images into the east-west direction is available (2002-08-19 06:10:41 UT).
In the text

	Fig. 3 Differential rotation profiles for the year 2014 (see Table 1). Sidereal rotation velocity is denoted by ω and heliographic latitude by b.
In the text

	Fig. 4 Differential rotation profiles calculated by different authors for the KSO data set and both hemispheres together (see corresponding rows in Table 3). Sidereal rotation velocity is denoted by ω and heliographic latitude by b.
In the text

	Fig. 5 Differential rotation profiles calculated by different authors for several data sets (GPR, EGR, Mt. Wilson, DPD, KSO) and both hemispheres together (see corresponding rows in Table 3). Sidereal rotation velocity is denoted by ω and heliographic latitude by b.
In the text

	Fig. 6 North – south asymmetry for the whole time period 1964–2016 and for DS and rLSQ methods (see also the corresponding rows in Table 2). Sidereal rotation velocity is denoted by ω and heliographic latitude by b. N – northern hemisphere, S – southern hemisphere.
In the text

	Fig. 7 Comparison between the present paper and Lustig (1983) differential rotation profiles calculated for solar cycle No. 20 (1964–1976) and for the northern and southern hemisphere separately, see corresponding rows in Tables 2 and 3. The sidereal rotation velocity is denoted by ω and heliographic latitude by b.
In the text